© ESO, 2012

1. Introduction

The interstellar medium (ISM) is in constant evolution, ruled by the energetic feedback from the cosmic cycle of star formation and stellar death. At the younger stages of star formation (bipolar outflows), and after the death of massive stars (SNRs), shock waves originating from the star interact with the ambient medium. They constitute an important mechanical energy input, and lead to the dispersion of molecular clouds and to the compression of cores, possibly triggering further star formation. Studying the signature of these interactions in the far-infrared and sub-mm range is paramount for understanding the physical and chemical conditions of the shocked regions and the large-scale roles of these feedback mechanisms.

Supernovae send shock waves through the ISM, where they successively carve out large hot and ionised cavities. They subsequently emit strong line radiations (optical/UV), and eventually interact with molecular clouds, driving lower-velocity shocks. Similar to their bipolar outflow equivalents, these shocks heat, compress, and accelerate the ambient medium before cooling down through molecular emission (van Dishoeck et al. 1993; Yuan & Neufeld 2011, hereafter Y11).

Valuable information has been provided by ISO (Cesarsky et al. 1999; Snell et al. 2005) and Spitzer (Neufeld et al. 2007, hereafter N07), but neither of those instruments provided sufficient spectral resolution to allow for a detailed study of the shock mechanisms. High-J CO emission is one of the most interesting diagnostics of SNRs. CO is indeed a stable and abundant molecule, and an important contributor to the cooling of these regions, whose high-frequency emission is expected to be a “pure” shock tracer. Observations of the latter must be carried out from above the Earth’s atmosphere. As part of a multi-wavelength study of MHD shocks that also includes Herschel data, we present here the first velocity-resolved CO (11–10) observations towards a prominent SNR-driven shock with the GREAT spectrometer onboard SOFIA, and combine them with new lower-J ones in a shock-model analysis.

2. The supernova remnant W28

W28 is an old (>10^4.5 yr, Claussen et al. 1999) SNR in its radiative phase of evolution, with a non-thermal radio shell centrally filled with thermal X-ray emission. Lying in a complex region of the Galactic disk at a distance of 1.9  ±  0.3 kpc (Velázquez et al. 2002), its structure in the 327 MHz radio continuum represents a bubble-like shape of about 40  ×  30 pc (Frail et al. 1993). Early on, molecular line emission peaks, not associated with star formation activity, but revealing broad lines, were suggested as evidence for interaction of the remnant with surrounding molecular clouds (Wootten 1981). Later studies spatially resolved the shocked CO gas layers from the ambient gas (Frail & Mitchell 1998; Arikawa et al. 1999). OH maser spots line up with the post-shock gas layers (Frail et al. 1994; Claussen et al. 1997; Hoffman et al. 2005), for which the strongest masers VLBA polarisation studies yield line-of-sight magnetic field strengths of up to 2 mG (Claussen et al. 1999). Pure rotational transitions of H₂ have been detected with ISO (Reach & Rho 2000) and were more recently observed with Spitzer, better resolved spatially and spectrally, by N07 and Y11.

Recently, very high energy (TeV) γ-ray emission has been detected by H.E.S.S. (Aharonian et al. 2008), Fermi (Abdo et al. 2010), and AGILE (Giuliani et al. 2010), spatially slightly extended and coincident with the bright interaction zones, W28-E and -F. If interpreted as the result of hadronic cosmic ray interactions in the dense gas (π⁰ decay), a cosmic ray density enhancement by an order of magnitude is required (which is supplied/accelerated by the SNR).

2.1. Sub-mm CO observations of W28F

APEX¹ (Güsten et al. 2006) observations towards W28F were conducted in 2009 and will be the subject of a forthcoming publication (Gusdorf et al., in prep.). For the present study, we used 100′′ × 100″ maps in the ¹³CO (3–2), CO (3–2), (4–3), (6–5), and (7–6) transitions, described in Appendix A.

Fig. 1 Overlay of the velocity-integrated CO (6–5) (colour background) with the CO (3–2) (white contours) emission observed with the APEX telescope. For both lines, the intensity was integrated between –30 and 40 km s-1. The wedge unit is K km s-1 in antenna temperature. The CO (3 − 2) contours are from 30 to 160 σ, in steps of 10σ = 16 K km s-1. The half-maximum contours of the CO (3 − 2) and (6–5) maps are indicated in red and black, respectively. The dark blue circle indicates the position and beam size of the SOFIA/GREAT observations. The APEX beam sizes of our CO (3 − 2), (4–3), (6–5), and (7 − 6) observations are also provided (upper right corner light green circles, see also Table A.1). The maps are centred at (RA[J2000] = , Dec[J2000]  =  −23°19′25″). The black and light blue hexagons mark the position of the OH masers observed by Claussen et al. (1997) and Hoffman et al. (2005).

In Fig. 1 the velocity-integrated CO (6–5) broad-line emission of W28F is shown overlaid with the CO (3–2) emission (white contours): a north-south elongated structure of about 100″ height and 30″ width traces the same warm accelerated post-shocked gas. In our high-resolution CO (6–5) data the structure is resolved, though probably still sub-structured similar to what is seen in H₂, e.g., Y11. Comparison with the distributions of excited H₂ and OH masers (whose locations also mark the leading edge of the non-thermal radio shell) suggests a textbook morphology of an SNR-molecular cloud interaction: the shock propagates E-NE into the ambient cloud that extends east for several arcmins. Hot H₂ and OH masers mark the first signposts of the shock-compressed gas. Farther downstream, the gas cooling is seen prominently in warm CO. The shock impact appears edge-on, but the fact that high – projected – streaming velocities are indeed observed (–30 km s^-1 with respect to the ambient cloud) requires a significant inclination angle.

Fig. 2 CO transitions observed in the position (+7″, –26″) indicated in Fig. 1: APEX (3–2), black (corresponding 13CO, green); (4–3), pink; (6–5), dark blue; (7–6), light blue; and SOFIA (11–10), red. The 13CO (3−2) and CO (11–10) profiles were multiplied by three for comparison purposes. Respective CO spectral resolutions are 0.212, 0.318, 0.159, 0.272, and 0.693 km s-1, and 0.664 km s-1 for 13CO (3–2).

We have selected the most prominent position in the southern extension of W28F, marked with the blue circle in Fig. 1, for our shock analysis. This position was also covered by Spitzer, offering a set of complementary H₂ data. Figure 2 shows the APEX spectra obtained towards this position, in main beam temperature units, all convolved to the 23.7″ beam of the SOFIA observations (Sect. 2.2). Absorption notwithstanding, the spectra show Gaussian profiles, with all line wings extending in blue (–30 km s^-1) and red (+15 km s^-1) velocities with respect to the ambient cloud (~10 km s^-1). The higher the line frequency, the narrower the line profile, with a typical linewidth of about 20 km s^-1 for CO (7–6). We detect red-shifted line-of-sight absorption features in all lines up to CO (6–5), the deepest arising at the ambient cloud velocity. Off-position contamination results in minor absorption features at 20 and 25 km s^-1 in the (3−2) and (4−3) profiles. Comparison with the ¹³CO(3−2) profile is also shown in Fig. 2. An analysis of the line temperature ratio of ¹²CO/¹³CO yields an optical thickness value of 3–7 in the wings of the ¹²CO (3–2), assuming a typical interstellar abundance ratio of 50–60 (e.g., Langer & Penzias 1993).

2.2. Far-infrared CO spectroscopy with GREAT/SOFIA

The observations towards W28F were conducted with the GREAT² spectrometer (Heyminck et al. 2012) during SOFIA’s flight from Stuttgart to Washington on September 21 2011. Only one position could be observed, towards the southern tip of the shocked cloud at offset (+7″, –26″) (Fig. 1). The CO (11–10) line was tuned to the frequency 1267.015 GHz LSB. The receiver was connected to a digital FFT spectrometer (Klein et al. 2012) providing a bandwidth of 1.5 GHz with a spectral resolution of 0.05 km s^-1. The observations were performed in double beam-switching mode, with an amplitude of 80″ (or a throw of 160″) at the position angle of 135° (NE-SW), and a phase time of 0.5 s. The nominal focus position was updated regularly against temperature drifts of the telescope structure. The pointing was established with the optical guide cameras to an accuracy of  ~5″. The beam width Θ_mb is 23.7″; the main beam and forward efficiencies are η_mb = 0.54 and η_f = 0.95. The integration time was 13 min ON source, for a final rms of 0.66 K. The data were calibrated with the KOSMA/GREAT calibrator (Guan et al. 2012), removing residual telluric lines, and subsequently processed with the CLASS software³.

The CO (11–10) spectrum, overlaid on the sub-mm lines in Fig. 2, reveals a markedly different profile: weak emission only is seen from the high-velocity gas, while the line is more prominent at low velocities. The profile basically follows the shape of the optically thin ¹³CO (3–2).

Fig. 3 Best-model comparisons between CO observations and models for each observed shock: the blue- (upper panel) and the red- (lower panel) components. Observations are marked by the black squares, our individual best-shock models are in blue line and circles (ζ1) or diamonds (ζ2), the warm layer that we used to compensate the ambient emission affecting the (3–2) and (4–3) transitions in green line and triangles, and the sum of each of these components is represented by a red line and circles or diamonds.

3. Discussion

3.1. The observations

Although most likely part of a single original shock clump, we separated the profiles into a blue lobe (–30 to 7.5–12.5 km s^-1) and a red lobe (7.5–12.5 to 40 km s^-1), and fitted the data independently. The uncertainty of the upper (lower) limit of those ranges reflects our lack of knowledge of the ambient velocity component. In Fig. 3 we plot the velocity-integrated intensities of the CO lines, all convolved to the same angular resolution of 23.7″, against the rotational quantum number of their upper level in a so-called “spectral line energy distribution” (or SLED, filled black squares with errorbars). The underlying assumption is that the filling factor is the same for all CO lines, which is validated by the similarity between their emitting regions (see for instance Fig. 1, where the half-maximum contours of CO (3–2) and (6–5) coincide at the available resolution). For Fig. 3, the assumption was made that all observed transitions present a circular emission region of radius 25″, corresponding to a filling factor of 0.53 (compatible with our maps, see Fig. 1).

The upper (lower) panel shows the diagram associated with the blue (red) shock component. The lower limit to the integrated intensity corresponds to the integration of the observed profile in the smallest velocity range, whereas upper limits were obtained by integrating Gaussian fits to the observation on the largest velocity range. The fits were adjusted to recover the shock flux lost through absorption and based on the un-absorbed parts of the profiles. Although yielding high errorbars on our measurements, specially for low J_up values, this method also provides the most conservative approach to our blue-red decomposition of the CO emission.

3.2. The models

We then compared the resulting integrated intensity diagrams to modelled ones. To build those, we used a radiative transfer module based on the large velocity gradient (LVG) approximation to characterise the emission from the CO molecule over outputs generated by a state-of-the-art model that calculates the structure of one-dimension, stationary shock layers (or approximations of non-stationary layers). This method has already been used and extensively introduced in Gusdorf et al. (2008). Since then, the LVG module has been modified to incorporate the latest collisional rate coefficients of CO with H₂ computed by Yang et al. (2010), but the shock model is the same. Based on a set of input parameters (pre-shock density, shock velocity, type, and age, and magnetic field parameter value b such as B[μG] = ), it calculates the structure of a shock layer, providing dynamical (velocity, density, temperature), and chemical (fractional abundances of more than 125 species linked by over 1000 reactions) variables values at each point of the layer. The relevant outputs are then used by our LVG module, which calculates the level population and line emissivities of the considered molecule. In the present case, the line temperatures were computed for transitions of CO up to (40–39) at each point of the shock layer, over which they were then integrated to form our model’s SLED. Because of the importance of H₂ (key gas coolant, abundant molecule), its radiative transfer was treated within the shock model. Its abundance was calculated at each point of the shock layer, based on the processes listed in Le Bourlot et al. (2002), and the populations of the 150 first levels were also calculated inside the shock code, their contribution to the cooling being included in the dynamical calculations. It is hence also possible to associate an H₂ excitation diagram to each model of the grid, see Sect. B.2.

Our model grid consists of a large sample of integrated intensity diagrams, obtained for stationary C- and J-type, and non-stationary CJ-type shocks. It covers four pre-shock densities (from 10³ to 10⁶ cm^-3), 0.45 ≤ b ≤ 2, over a range of velocities spanning from 5 km s^-1 up to the maximum, “critical” value above which a C-shock can no longer be maintained, which depends on the other parameters. We then independently compared the observed blue- and the red-shock component to the whole grid diagrams. We used a χ² routine that was set up to provide the best fits to our purest shock-tracing lines: CO (6–5), (7–6), and (11–10), to reduce any effect induced by ambient emission contamination on the lower-lying transitions.

3.3. The results

The upper (lower) panel of Fig. 3 shows our best-fitting models for the blue (red) component. We found that only stationary, C-type shock models can reproduce the levels of observed CO (6−5), (7–6) and (11–10) integrated intensities. In both cases, the pre-shock density is 10⁴ cm^-3, for which a stationary state is typically reached for a shock age of 10⁴ years, in agreement with the age of the remnant of 3.5 × 10⁴ years quoted by Giuliani et al. (2010). The modelled shock velocities, 25 and 20 km s^-1, respectively, are rough upper limits to the observed ones, an expected conclusion given the one-dimensional nature of our models. Additionally, our models provide respective constraints to the magnetic field component perpendicular to the shock layers, with parameter b of 1 and 0.45. The difference between those two values can be explained by a projection effect. In both cases, the direct post-shock density is of the order of 2–3 × 10⁵ cm^-3, yielding expected post-shock values for this component of 0.05–0.55 mG, compatible with those inferred by Claussen et al. (1999) (2 mG for the total post-shock magnetic field), and by Hoffman et al. (2005) (0.3–1.1 mG based on the OH maser measurements in the positions indicated in Fig. 1 located within the blue beam). Finally, our models produce typical OH column densities of a few 10¹⁶ cm^-2. This complies with the list of requirements set out by Lockett et al. (1999) for exciting the observed OH 1720 MHz masers, together with the C-nature of the considered shocks, as well as the temperature and density conditions they create.

The respective CO column density found in each of the red and blue shock models is 11 and 5 × 10¹⁷ cm^-3. However, our shock models failed to account for the important observed levels of CO (3–2) and (4–3) emission. We found it possible to model this lower-J surplus emission by adding a thin layer of gas to our shock model results. Its corresponding emission was calculated with our LVG module used in a “homogeneous slab” mode, and indeed generated significant emission only from the (3–2) and (4–3) lines, as can also be seen in Fig. 3. We used a linewidth of 10 km s^-1 and a density of n_H = 10⁴ cm^-3 with a fractional abundance of X(CO) = 10^-4 for the blue and red cases. We found that this “CO layer” is warm (35 K and 15 K for the blue and red components, respectively), and adds to the pure-shock CO column density in both components, from 2 × 10¹⁷ cm^-2 for the blue component, and 10¹⁷ cm^-2 for the red one. Our final result, the sum of the C-shock and “CO layer”, is also shown in Fig. 3.

To understand the physical origin of this emitting layer, we varied the cosmic ray ionisation rate in our best-fitting models from its solar value ζ₁ = 5 × 10^-17 s^-1 to the ζ₂ = 3.4 × 10^-16 s^-1 indicated by Hewitt et al. (2009) in the specific case of W28. In our modelling, this modification impacts on the chemistry (higher ionisation fraction) and the physics (warmer gas) of our outputs. The resulting shock profile is more narrow (hence generating slightly less CO emission, see Fig. 3), with higher post-shock temperatures close to the aforementioned CO layer.

Because one of the limitations of our shock model is the absence of significant velocity gradient in the post-shock region, the warm CO layer component might correspond to a post-shock layer of gas mildly heated by the energetic radiation that seems to exist in the region (UV or more energetic, as characterised by Hewitt et al. 2009). The rigorous inclusion of these radiative effects in our shock model is still in progress. In Sect. B we present a consistency check based on the use of the Spitzer H₂ observations.

Online material

Appendix A: The APEX observations

Fig. A.1 Location of the field covered by our CO observations on the larger-scale radio continuum image at 327 MHz taken from Claussen et al. (1997): entire SNR in the left panel, zoom in the right panel. The shown CO observations are the (6–5) map, also displayed in Fig. 1.

APEX observations towards the supernova remnant W28F were conducted in several runs in the year 2009 (in May, June, August and October). We used of a great part of the suite of heterodyne receivers available for this facility: APEX–2 Risacher et al. (2006), FLASH460 Heyminck et al. (2006), and CHAMP⁺Kasemann et al. (2006); Güsten et al. (2008), in combination with the MPIfR fast Fourier transform spectrometer backend (FFTS, Klein et al. 2006). The central position of all observations was set to be , β_[J2000]  =  −23°18′5850. Focus was checked at the beginning of each observing session, after sunrise and/or sunset on Mars, or on Jupiter. Line and continuum pointing was locally checked on RAFGL1922, G10.47B1, NGC 6334-I or SgrB2(N). The pointing accuracy was found to be of the order of 5″ rms, regardless of the receiver that was used. Table A.1 contains the main characteristics of the observed lines and corresponding observing set-ups: frequency, beam size, sampling, used receiver, observing days, forward and beam efficiency, system temperature, spectral resolution, and finally the velocity interval that was used to generate the integrated intensity maps. The observations were performed in position-switching/raster mode using the APECS software Muders et al. (2006). The data were reduced with the CLASS software (see http://www.iram.fr/IRAMFR/GILDAS). For all observations, the maximum number of channels available in the backend was used (8192), except for CO (4–3), for which only 2048 channels were used, leading to the spectral resolutions indicated in Table A.1. Maps were obtained for all considered transitions, covering the field introduced in Fig. 1 and put in the perspective of the whole SNR in Fig. A.1.

Appendix B: The H₂ observations

B.1. The dataset

Fig. B.1 Overlay of the map of CO (6–5) emission observed by the APEX telescope (colour background) with the H2 0–0 S(5) emission (white contours), observed with the Spitzer telescope. The wedge unit is K km s-1 (antenna temperature) and refers to the CO observations. The H2 0–0 S(5) contours are from 50 to 210σ, in steps of 20σ ≃ 1.6 × 10-4 erg cm-2 s-1 sr-1. The green contour defines the half-maximum contour of this transition. Like in Fig. 1, the blue circle indicates the beam size of the SOFIA/GREAT observations, on the position of the centre of the circle, marked by a blue dot. The beam and pixel sizes of the CO (6–5), and H2 0–0 S(5) observations are also provided in yellow (big and small, respectively) circles in the lower right corner. The black contour delineates the half-maximum contour of the H2 0–0 S(2) transition. The field is smaller than in Fig. 1, and the (0, 0) position is that of the Spitzer observations (see Appendix B).

As a consistency check for our models, we used of the Spitzer/IRS observations of the H₂ pure rotational transitions (0–0 S(0) up to S(7)), reported and analysed in N07 and Y11. Although the original dataset also includes other ionised species, we chose to use only the H₂ data. The raw product communicated to us by David Neufeld contains rotational transitions maps, with 1.2″ per pixel, centred on , β_[J2000]  =  −23°19′2492. Figure B.1 shows an overlay of our APEX CO (6–5) map with the H₂ 0–0 S(5) region observed by Spitzer. The figure shows coinciding maxima between the two datasets in the selected position, and a slightly different emission distribution. This might be the effect of the better spatial resolution of the H₂ data, which reveal more peaks than in CO. This overlay also shows the slightly different morphology of the emission of the S(2) (half maximum contour in black) and S(5) (half maximum contour in green) transitions, at the available resolution.

B.2. Excitation diagram

We performed a consistency check on our modelling (see Sect. 3) using the excitation diagram derived for the selected emission region. The H₂ excitation diagram displays ln(N_νj/g_j) as a function of E_νj/k_B, where N_νj (cm^-2) is the column density of the rovibrational level (v,J), E_νj/k_B is its excitation energy (in K), and g_j = (2j + 1)(2I + 1) its statistical weight (with I = 1 and I = 0 in the respective cases of ortho- and para-H₂). If the gas is thermalised at a single temperature, all points in the diagram fall on a straight line. The selected position for the present study was introduced in Sect. 2.1, and can be seen in Figs. 1 and B.1. The column density of the higher level of each considered transition was extracted by averaging the line intensities in a 11.85″ radius circular region, consistent with our handling of the CO data. In the process, we corrected the line intensities for extinction, adopting the visual extinction values from N07, A_ν = 3 − 4, and using the interstellar extinction law of Rieke & Lebofsky (1985). The resulting excitation diagram can be seen in Fig. B.2. We note that the intensity of the H₂ 0–0 S(6) transition (J = 8) cannot be determined reliably, because the line is blended with a strong 6.2 μm PAH feature.

B.3. Comparisons with our models

Fig. B.2 H2 excitation diagram comparisons. Upper panel: evolution of the modelled excitation diagram obtained for ζ1 = 5 × 10-17 s-1, varying the initial value of the OPR from the unrealistic, extreme-case 0 value to its equilibrium one, 3. Lower panel: influence of the cosmic ray ionisation rate variation on our modelled excitation diagram, with the initial OPR set to 3.

Unlike our CO observations, the Spitzer/IRS observations are not spectrally resolved (owing to the relatively low resolving power of the short low module, in the range 60–130). On the other hand, given the minimum energy of the levels excited in the H₂ transitions (for 0–0 S(0), E_u ~ 509.9 K), we can make the double assumption that the lines are not contaminated by ambient emission or self-absorption, and that the measured line intensities are the result from the blue- and red-shock present in the line of sight. Therefore, one must compare the observed H₂ level populations to what is generated by the sum of our best-fitting blue- and red-shock models. Figure B.2 shows the comparison between the excitation diagram derived from the sum of our two CO best-fitting models and the observed one. The H₂ excitation diagram comprising the sum of our CO best-fitting models is shown in both panels in red circles, and provides a satisfying fit to the observations. With the aim of improving the quality of this fit, we also studied the influence of the initial ortho-to-para ratio (OPR) value, upper panel, and of the cosmic ray ionisation rate value, lower panel.

In the upper panel, for a cosmic ray ionisation rate value adopted as the solar one, ζ₁ = 5 × 10^-17 s^-1, the influence of the initial OPR value is studied. In their models, N07 inferred a value of 0.93 for the “warm” component (~322 K). On the other hand, the value associated to their “hot” (~1040 K) component could not be estimated, owing to the large uncertainty associated to the determination of the H₂ 0–0 S(6) flux. In our shock models, the OPR value is consistently calculated at each point of the shocked layer, and is mostly the result of conversion reactions between H₂ and H, H⁺ or H. The heating associated to the passage of the wave increases the OPR value towards the equilibrium value of 3.0. Nevertheless, this value is reached only when it is also the initial one. The N07 situation is then probably adequately approximated in our models where the OPR initial value is less than 3. However, the influence of this parameter is minimal on the excitation diagram, where it only seems to create a slight saw-tooth pattern between the odd and even values of J.

We finally investigated the effects of the high-energy radiation field on the excitation diagram in the lower panel of Fig. B.2. The limits of our modelling are reached because the only handle that we have to study this effect on our outputs is the variation of the cosmic ray ionisation rate value, which affects the corresponding chemistry (higher ionisation fraction) and physics (warmer gas). A proper treatment of the energetic radiation components listed in Hewitt et al. (2009) should indeed incorporate UV pumping or self-shielding, as would be the case in PDR regions (e.g., Habart et al. 2011), as well as chemical and H₂ excitation effects by X-ray (see Dalgarno et al. 1999) and cosmic rays (e.g., Ferland et al. 2008). Keeping these limits in mind, we found that even strong modifications of the cosmic ray ionisation rate from the solar value of 5 × 10^-17 s^-1 to the more extreme value of 10^-14 s^-1 have only very limited effects on the excitation diagram. The rather convincing fits to the H₂ data seem to indicate that these excitation effects by energetic photons could well be minimal in our case, but we repeat that their proper inclusion to our models is work in progress.

Acknowledgments

We thank the SOFIA engineering and operations teams, whose support has been essential during basic science flights, and the DSI telescope engineering team. Based [in part] on observations made with the NASA/DLR Stratospheric Observatory for Infrared Astronomy. SOFIA Science Mission Operations are conducted jointly by the Universities Space Research Association, Inc., under NASA contract NAS2-97001, and the Deutsches SOFIA Institut under DLR contract 50 OK 0901. S. Anderl acknowledges support by the Deutsche Forschungsgemeinschaft within the SFB 956 “Conditions and impact of star formation”, the International Max Planck Research School (IMPRS) for Radio and Infrared Astronomy at the Universities of Bonn and Cologne, and the Bonn-Cologne Graduate School of Physics and Astronomy. A. Gusdorf acknowledges support by the grant ANR-09-BLAN-0231-01 from the French Agence Nationale de la Recherche as part of the SCHISM project.

References

All Tables

All Figures

Fig. 1 Overlay of the velocity-integrated CO (6–5) (colour background) with the CO (3–2) (white contours) emission observed with the APEX telescope. For both lines, the intensity was integrated between –30 and 40 km s-1. The wedge unit is K km s-1 in antenna temperature. The CO (3 − 2) contours are from 30 to 160 σ, in steps of 10σ = 16 K km s-1. The half-maximum contours of the CO (3 − 2) and (6–5) maps are indicated in red and black, respectively. The dark blue circle indicates the position and beam size of the SOFIA/GREAT observations. The APEX beam sizes of our CO (3 − 2), (4–3), (6–5), and (7 − 6) observations are also provided (upper right corner light green circles, see also Table A.1). The maps are centred at (RA[J2000] = , Dec[J2000]  =  −23°19′25″). The black and light blue hexagons mark the position of the OH masers observed by Claussen et al. (1997) and Hoffman et al. (2005).

In the text

Fig. 2 CO transitions observed in the position (+7″, –26″) indicated in Fig. 1: APEX (3–2), black (corresponding 13CO, green); (4–3), pink; (6–5), dark blue; (7–6), light blue; and SOFIA (11–10), red. The 13CO (3−2) and CO (11–10) profiles were multiplied by three for comparison purposes. Respective CO spectral resolutions are 0.212, 0.318, 0.159, 0.272, and 0.693 km s-1, and 0.664 km s-1 for 13CO (3–2).

In the text

Fig. 3 Best-model comparisons between CO observations and models for each observed shock: the blue- (upper panel) and the red- (lower panel) components. Observations are marked by the black squares, our individual best-shock models are in blue line and circles (ζ1) or diamonds (ζ2), the warm layer that we used to compensate the ambient emission affecting the (3–2) and (4–3) transitions in green line and triangles, and the sum of each of these components is represented by a red line and circles or diamonds.

In the text

	Fig. A.1 Location of the field covered by our CO observations on the larger-scale radio continuum image at 327 MHz taken from Claussen et al. (1997): entire SNR in the left panel, zoom in the right panel. The shown CO observations are the (6–5) map, also displayed in Fig. 1.
In the text

Fig. B.1 Overlay of the map of CO (6–5) emission observed by the APEX telescope (colour background) with the H2 0–0 S(5) emission (white contours), observed with the Spitzer telescope. The wedge unit is K km s-1 (antenna temperature) and refers to the CO observations. The H2 0–0 S(5) contours are from 50 to 210σ, in steps of 20σ ≃ 1.6 × 10-4 erg cm-2 s-1 sr-1. The green contour defines the half-maximum contour of this transition. Like in Fig. 1, the blue circle indicates the beam size of the SOFIA/GREAT observations, on the position of the centre of the circle, marked by a blue dot. The beam and pixel sizes of the CO (6–5), and H2 0–0 S(5) observations are also provided in yellow (big and small, respectively) circles in the lower right corner. The black contour delineates the half-maximum contour of the H2 0–0 S(2) transition. The field is smaller than in Fig. 1, and the (0, 0) position is that of the Spitzer observations (see Appendix B).

In the text

	Fig. B.2 H2 excitation diagram comparisons. Upper panel: evolution of the modelled excitation diagram obtained for ζ1 = 5 × 10-17 s-1, varying the initial value of the OPR from the unrealistic, extreme-case 0 value to its equilibrium one, 3. Lower panel: influence of the cosmic ray ionisation rate variation on our modelled excitation diagram, with the initial OPR set to 3.
In the text