Constraining physical conditions for the PDR of Trumpler 14 in the Carina Nebula

We investigate the physical conditions of the CO gas near the young star cluster, Trumpler 14 of the Carina Nebula. The observations presented in this work are taken with the Fourier Transform Spectrometer (FTS) of the Spectral and Photometric Imaging REceiver (SPIRE) onboard the Herschel Space Observatory. Our field of view covers the edge of a cavity carved by Trumpler 14 about $1\,\mathrm{Myr}$ ago and marks the transition from HII regions to photo-dissociation regions. With the state-of-the-art Meudon PDR code, we successfully derive the physical conditions, which include the thermal pressure ($P$) and the scaling factor of radiation fields ($G_{\mathrm{UV}}$), from the observed CO spectral line energy distributions~(SLEDs) in the observed region. The derived $G_{\mathrm{UV}}$ values generally show an excellent agreement with the UV radiation fields created by nearby OB-stars and thus confirm that the main excitation source of the observed CO emission are the UV-photons provided by the massive stars. The derived thermal pressure is between $0.5-3\,\times\,10^{8}\,\mathrm{K\,cm^{-3}}$ with the highest values found along the ionization front in Car I-E region facing Trumpler 14, hinting that the cloud structure is similar to the recent observations of the Orion Bar. Comparing the derived thermal pressure with the radiation fields, we report the first observationally-derived and spatially-resolved $P \sim 2\times10^4\,G_{\mathrm{UV}}$ relationship. As direct comparisons of the modeling results to the observed $^{13}\mathrm{CO}$, [OI] $63\,\mathrm{\mu m}$, and [CII] $158\,\mathrm{\mu m}$ intensities are not straightforward, we urge the readers to be cautious when constraining the physical conditions of PDRs with combinations of $^{12}\mathrm{CO}$, $^{13}\mathrm{CO}$, [CI], [OI] $63\,\mathrm{\mu m}$, and [CII] $158\,\mathrm{\mu m}$ observations.


Introduction
The radiative feedback of stars on their parent cloud is a key topic both in the context of evolution and composition of interstellar matter and to constrain star formation mechanisms. Emission lines from photo-dissociation regions (PDRs) have been studied for decades to understand the physical and chemical processes induced by this feedback (Tielens & Hollenbach 1985;Sternberg & Dalgarno 1989). In these regions, far ultraviolet (FUV) photons dissociate molecules and heat the gas at several hundred Kelvins via photo-electric effect on grains and, in dense gas, by H 2 pumping followed by collisional de-excitations (Röllig et al. 2006). Several tracers are commonly used to study PDRs. Observations of atomic lines of O, C + and C probe the neutral hydrogen layer and give access to most of the gas cooling rate at the edge of the cloud. The atomic-molecular transition can be observed thanks to the emission of H 2 infrared lines. Less abundant molecules, diatomic and complex ones such as carbon chains, formaldehyde or methanol, have been detected at various depths in PDRs showing that complex chemistry takes place in this hostile environment (Young Owl et al. 2000;Lis & Schilke 2003;Guzmán et al. 2013Guzmán et al. , 2014Nagy et al. 2013;Cuadrado et al. 2015). Because H 2 is difficult to observe, except at the warm surface layer of cloud due to its lack of electric dipole moment, many studies rely on the observation of CO and its A&A 618, A53 (2018) isotopologs to constrain the mass of gas and the morphology of clouds.
Since its launch in 2009, Herschel Space Observatory (Pilbratt et al. 2010) has allowed Galactic and extragalactic observations of excited CO in high-J states (J up > 10) in PDRs of various objects (Köhler et al. 2014;Stock et al. 2015;Wu et al. 2015;Parikka et al. 2017). In parallel, observations of CO emission in protostars and extragalactic environments showed highly excited CO. For example, emission from rotational levels J up > 40 has been observed in the protostellar region Orion BL/KL, suggesting an excitation by outflows and shocks (Goicoechea et al. 2015). In extragalactic regions, the origin of CO excitation is unclear (Yıldız et al. 2010;Hailey-Dunsheath et al. 2012;Greve et al. 2014;Kamenetzky et al. 2016;Rosenberg et al. 2015). Candidate excitation sources include mechanical-heating, shocks, X-rays, and UV photons, depending on the objects and the adopted tools in the analyses. Whether it is possible to infer the UV radiation field intensity and star formation rate from high-J CO line intensities has been a topic of interest. The study of CO excitation in spatially resolved Galactic PDRs can help to answer this question. Presently, most analyses rely on the usage of non-local thermal equilibrium (non-LTE) radiative transfer codes, such as large velocity gradient (LVG) approximations, to infer gas density and temperature from CO line intensities (Köhler et al. 2014). Very few detailed analyses based on PDR modeling, which constrains the thermal and chemical processes hidden behind the observed CO lines, can be found in the literature. Based on the non-spatially resolved CO observations from J = 4−3 up to J = 19−18 in NGC 7023 and up to J = 23−22 in the Orion Bar, as well as other atomic and molecular lines, and their comparisons with the Meudon PDR code (Le Petit et al. 2006), it is suggested that the FUV radiation field produced by the Trapezium star cluster and HD 200775 are responsible for a high pressure layer (P ∼ 10 8 K cm −3 ) at the very edge of the two PDRs and is sufficient to support the high-J CO excitation observed in the two objects (Joblin et al. 2018). These results are in agreement with the recent ALMA observations which confirm the presence of a thin layer with high pressure at the edge of the PDR in the Orion Bar (Goicoechea et al. 2016. Further such studies with detailed PDR modeling are necessary for understanding the excitation of CO in the interstellar environment. Another well-known object that is suitable for studying the impact of radiative feedback of massive stars on their parent cloud is the Carina Nebula. Located at a close distance, 2.3 kpc (Smith 2006b), from us in the Milky Way, the Carina Nebula is the largest and brightest nebula in the southern sky. Some of its central regions are even brighter than the famous Orion Nebula. Aside from the famous luminous blue variable, η Carinae (η Car), the Carina Nebula hosts 70 known O stars, 3 Wolf-Rayet stars, and 127 B0 to B3 stars (Gagné et al. 2011).
It was first revealed in the observations of radio continuum at 6 cm that the HII regions of the Carina Nebula include two major components, Car I and Car II, separated by ∼10 (∼7 pc) along the Galactic plane (McGee & Gardner 1968;Gardner & Morimoto 1968). Based on the distribution of local maxima of radio continuum brightness at 0.843 GHz, Whiteoak (1994) detailed the structure of Car I into Car I-W, Car I-E, and Car I-S, which are associated with the ionization front shaped by the UV-photons from Trumpler 14 in the west, east and south of Car I. With the NASA Kuiper Airborne Observatory (KAO) far-infrared (FIR) and Fleurs Synthesis Telescope (FST) 1415 MHz observations, it is confirmed that the main heating sources for Car I and Car II are Trumpler 14 and Trumpler 16, respectively (Harvey et al. 1979;Retallack 1983). Figure 1 gives an overview of the north region of the Carina Nebula observed by the UK Schmidt Telescope in the R-band of the Second Digital Sky Survey (DSS2 red) 1 .
In this paper, we study the CO and C emission observed by Herschel SPIRE/FTS across an HII-PDR interface located at the northwest in the Carina Nebula to better understand the origin of CO excitation in PDRs as well as the impact of radiative feedback of massive stars. The observed region (enclosed by the white square in Fig. 1) covers Car I-E (upper right), Car I-S (center), and a small region at the intersection of Car I and Car II (Car I/II, lower left, see also Fig. 2 for the region definition). Our field-of-view (FOV) covers the edge of a cavity likely shaped by the young (1−2 Myr) OB-star cluster Trumpler 14 since ∼1 Myr ago with its ionization front sitting at a projected distance of ∼2 pc from Trumpler 14.
Trumpler 14 hosts a dozen O-stars and more than a hundred B-stars. Its stellar components have been well studied (Vazquez et al. 1996;Carraro et al. 2004). The estimated UV input by the stellar members of Trumpler 14 to the ionization front of PDR is 10 4 G 0 (Brooks et al. 2003) where G 0 = 1.6 × 10 −3 erg s −1 cm −2 according to Habing (1968, hereafter, the Habing unit). For comparison, this FUV flux is slightly lower than the one induced by the Trapezium star cluster on the famous Orion Bar (∼3 × 10 4 G 0 , Goicoechea et al. 2015) and comparable to the one found in extragalactic starburst regions (Malhotra et al. 2001;Chevance et al. 2016). Even though the FUV flux is similar between the Carina Nebula and the Orion Bar, the distance of the Trapezium star cluster to the Bar is approximately ten times closer than in the Carina Nebula. As Trumpler 14 is a more evolved star cluster than the Trapezium star cluster, the interaction zones at its PDR interface can be better approximated by pressure equilibrium in a stationary PDR model (Brooks et al. 2003). Within the dense molecular gas at the northwest region, young stellar objects (YSOs) are carving out ionized regions anew, proving the on-going and active star formation in this region (Tapia et al. 2015). Such a region is likely a representation of prototypical PDRs near star-forming sites that dominate the extragalactic observations of starburst systems where individual PDRs are mostly unresolved in the infrared and submm.
The PDRs in the Car I-E region have thus far been studied by several groups. The FUV flux in this region, as estimated from stellar composition and FIR observations, is ∼10 4 in the Habing unit (Brooks et al. 2003;Mizutani et al. 2004). However, estimations with PDR models have led to lower values of FUV flux (1390, Oberst et al. 2011 and3200, Kramer et al. 2008 in the Habing unit). Looking at the bigger picture, these results indicate that the UV-photons provided by nearby massive stars are more than sufficient in supporting the global heating of PDRs observed in this region. However, the inconsistency persists between the model-estimated values of FUV flux and  (Tapia et al. 2003). The areas enclosed by the white dotted and dashed lines indicate the regions previously observed by Brooks et al. (2003) and Kramer et al. (2008). The long dashed line across the image marks the Galactic plane at b = −0.72 • . The square enclosed by the solid white line indicates the FOV in this work. It spans an area of 7 × 7 (4.7 × 4.7 pc). A zoom-in view of the area included in the white square can be found in Fig. 2. the UV radiation fields estimated from stellar composition and FIR observations.
In this work, we adopt a distance of 2.35 ± 0.5 kpc to the Carina Nebula (Smith 2006b), although it is worth noting that the photometric analysis from pre-main sequence (PMS) stars in the nebula gives a distance of ∼2.9 kpc to the star clusters, Trumpler 14 and 16 (Hur et al. 2012). All the sky coordinates quoted in this work are in the J2000 epoch. In the model analysis, we adopt the scaling factor according to the radiation field given by Mathis et al. (1983, hereafter, the Mathis unit) unless otherwise specified. The difference between the Habing and Mathis units is about a factor of 1.3, in the spectral range of 91.2 to 111 nm. Our observations and data reduction procedure are presented in Sect. 2. We interpret these data with the Meudon PDR code in Sect. 3. Results are then presented and discussed in Sects. 4 and 5, respectively.

Observations
The Herschel data presented in this work have been observed on 2012-12-03 and 2013-02-04 as one of the Open Time (OT1) programs (PI: Takashi Onaka). Our results are mainly derived from the spectroscopic data observed by the Fourier Transform Spectrometer (FTS) of the Spectral and Photometric Imaging REceiver (SPIRE, Griffin et al. 2010) on board the Herschel Space Observatory (Pilbratt et al. 2010). We give a brief instrumentation overview of the Herschel SPIRE/FTS, and then layout the steps we take to map the observations. A&A 618, A53 (2018) Fig. 2. A zoom-in of our FOV (the white square in Fig. 1). The effective FOV is indicated by the pixels in gray-colors. The background image is the optical image taken in the DSS2-red band. The positions of the central bolometer, SLWC4, at the six requested pointings are indicated with the magenta crosses and labeled accordingly as in Table 1. The boundaries of Car I-E, Car I-S, and Car I/II areas are indicated with the thick white lines. The locations of two massive stars, ALS 15 204 and ALS 15 203, discussed in Sect. 4.2 are marked with the cyan and yellow "×" in the figure. The CO SLED observed from the three yellow-masked pixels are displayed in Fig. 8 and discussed in Sect. 5.1.

Overview of the Herschel SPIRE/FTS
The Herschel SPIRE/FTS is made of two bolometer arrays, the SPIRE Short Wavelength (SSW, and SPIRE Long Wavelength (SLW, 316−672 µm) arrays. The SSW and SLW arrays pack 19 and 7 unvignetted bolometers, respectively, hexagonally onto the focal plane and each simultaneously cover a 2.6 diameter of field of view (FOV) on the sky. Our observations were completed in six single pointings on the operational days 1299 and 1362. At each single-pointing position, we observe the sky in the intermediate spatial-sampling (4 jiggle pointings) and high-spectral-resolution (∆σ = 0.04 cm −1 or ∆ν = 1.2 GHz) modes. The information of our observation is listed in Table 1. For each jiggle pointing, the sky was observed in eight scans, which add up to an effective exposure time of approximately 266 s across the FOV. The spectra presented in this work are calibrated against Uranus, i.e. point-source calibrated, using the SPIRE calibration tree, spire_cal_14_3. Details of the SPIRE/FTS calibration can be found in Swinyard et al. (2014). Beginning with the unapodized spectra, averaged from all scans, and their standard deviation used as the uncertainties, we describe our map-making strategy in the following subsection.

Making maps of the Herschel SPIRE/FTS lines
In the intermediate-sampling mode, the bolometer spacings for the SSW and SLW arrays are 16.3 and 25.3 , respectively, which are approximately the sizes of individual beams at their minimum within the bandwidths. Although the individual single-pointing observation used in this work is not fully Nyquist sampled, combining the footprints of bolometers from all six  Fig. 2. The Obs_ID is the unique identifier for Herschel observations. The labels used in Onaka et al. (2008) are given in the Note column.  Table 2, which correspond to 8.3 and 16.5 for the SSW and SLW arrays, respectively.
single-pointing positions (see Fig. 2), the average bolometer spacings in our FOV achieve 8.6 and 12.2 for SSW and SLW arrays, respectively. Histograms of distances between any given bolometer to its closest neighbor for SSW and SLW arrays are shown in Fig. 3, and they demonstrate that the observation within the effective FOV is approximately fully Nyquist sampled.
Within the bandwidth (194 < λ < 672 µm) of the Herschel SPIRE/FTS, the relevant atomic and molecular lines used in this study are listed in Table 2. As discussed in Wu et al. (2013) and Swinyard et al. (2014), the beam full-width-halfmaxima (FWHM) and sensitivities of bolometers vary across the instrument bandwidth. At each observed emission line, we list the corresponding instrument FWHM and mean 1σ uncertainty of the integrated line intensity ( ) in Table 2.
Around the frequency (ν line − 15 GHz < ν < ν line + 15 GHz) of each line listed in Table 2, the observed spectra from all unvignetted bolometers are assigned to a regular grid map. The world coordinate system (WCS) of the maps is centered at (α,δ) = (10:43:45.07, -59:36:19.63) with each pixel sized 15 × 15 (0.17 × 0.17 pc) and a total dimension of 7 × 7 (4.7 × 4.7 pc). We follow a similar map-making recipe to that detailed in Sect. 2.1.2 of Wu et al. (2015) and briefly discuss the main steps here. The stand-alone map-making software, HFTS_mapping, is written in Python and uses utilities publicly available in astropy 2 (version 1.1.1, Robitaille et al. 2013), scipy (version 0.17.0) and numpy (version 1.12.1) 3 .
HFTS_mapping is now available for public use 4 .
Our map-making procedure can be laid out in three steps: -First, before coadding the level-two bolometer spectra on the map, a continuum level is estimated between two frequency intervals, (ν line − 15 GHz < ν < ν line − 5 GHz) and (ν line + 5 GHz < ν < ν line + 15 GHz), and removed as a parabola function. Weighted by the uncertainties, which are estimated as the standard deviation from all scans, the continuum-removed bolometer spectra are then assigned to the pixels, according to the bolometer pointings, to create spectral cubes over the FOV for each line listed in Table 2. -Second, we measure the integrated line intensity, assuming a sinc profile, from spectral cubes constructed in the previous step. For CO J = 7−6 and [CI] 3 P 1 − 3 P 2 , which are only 2.7 GHz apart, their integrated line-intensities are measured simultaneously with two sinc profiles whose central frequencies are fixed to differ by 2.7 GHz. At the end of this step, the maps of integrated line-intensities are produced (see figures in Appendix A). -The last step of map-making is to convolve the integrated line-intensity maps to the same spatial resolution, 42 (∼0.5 pc), which is the largest beam FWHM among all observed lines (see Table 2). We refer the curious readers to Wu et al. (2015) for a more detailed description on the generation of kernels used in convolution. Uncertainties for the final integrated line-intensities are estimated by a Monte-Carlo experiment which comprises 300 repetitions of the above three steps. For each repetition, the uncertainties used for coaddition of bolometer spectra in the first step are replaced by randomly generated noises which include random and systematic uncertainties from the observation. The random uncertainties are represented by a normal distribution whose standard deviation is the root-mean-square (rms) of residuals from a given bolometer spectrum after subtracting its continuum and the sinc line-profile(s).
For all bolometer spectra taken in the same single-pointing observation, one normal distribution that represents the systematic uncertainty has also been taken into account. The 1-σ of the normal distribution for the systematic uncertainty is set to be 10% of the bolometer spectrum (Swinyard et al. 2014). The uncertainties used in the first step are then replaced by the rms of A&A 618, A53 (2018) the generated random and systematic uncertainties. At the end of the Monte Carlo experiment, 300 convolved maps of integrated line-intensities for each line are produced. For a given line on each pixel, the standard deviation between the 300 maps is taken as the final uncertainty. Within our FOV, the CO J = 4−3 and J = 7−6 transitions have previously been observed with the 4 m NANTEN2 telescope. The FWHMs of the NANTEN2 observations are 38 and 26.5 for the CO J = 4−3 and J = 7−6 transitions, respectively, with 20% calibration uncertainty (Kramer et al. 2008). After convolving the CO J = 4−3 and J = 7−6 images observed by NANTEN2 to the instrument FWHM of the Herschel SPIRE/FTS, the data observed by the two telescopes show agreement within 30%, while the observed values by the Herschel SPIRE/FTS are systematically higher. Taking the uncertainty in the atmospheric transmission model (Pardo et al. 2001) into account, the observations compare reasonably well.

Data interpretation
We deduce the physical conditions at each pixel with version 1.5.2 of the Meudon PDR code (Le Petit et al. 2006), which is available on the ISMServices webpage 5 of the Paris Observatory. In this section, we first give an overview of the Meudon PDR code and discuss the strategies adopted in the data interpretation.

The Meudon PDR code
The Meudon PDR code computes chemical densities and gas excitation in a one-dimensional (1D) stationary plane-parallel slab of gas and dust illuminated by a radiation field. At each point of the spatial grid, the code consistently solves the radiative transfer, the chemistry and the thermal balance, which computes the heating and cooling mechanisms individually. The chemical network used in this paper includes 224 species and 6311 chemical reactions.
H 2 formation on grains is a key process since it controls the transition between the neutral atomic and the molecular phases. The default prescription in the Meudon PDR code simulates the H 2 formation on dust grains through Eley-Rideal and Langmuir-Hinshelwood mechanisms, as described in Le Bourlot et al. (2012). A more complex yet computing-time-demanding H 2 formation prescription that simulates the stochastic heating of grains by UV photons and its effect on the H 2 formation and ortho/para conversion is also available for the Meudon PDR code (Bron et al. 2014(Bron et al. , 2016. In order to gain efficiency in computing time, we opt for the Le Bourlot et al. (2012) prescription.
At each point of the spatial grid, the radiative transfer equation is solved considering absorption and diffusion by dust as well as absorption in continuum by some gas species, such as carbon. Self-shielding for critical species, such as H 2 and CO, is computed with the Federman et al. (1979) formalism (see Le Petit et al. 2006 for details).
The most important heating processes taken into account are the photo-electric effect on grains (Bakes & Tielens 1994), heating by collisional de-excitation of H 2 , and cosmic rays as well as exothermic reactions. As for the cooling processes, the non-LTE level excitation of main coolants is computed at each position in the cloud, considering various micro-physical processes, such as collisions and radiative (de-)excitation (Gonzalez Garcia et al. 2008).
Another important aspect is the treatment of grains. Grains impact the photo-electric heating rate, the penetration of UV photons, the IR dust emission, and the formation of H 2 . In the models used for this work, a Mathis-Rumpl-Nordsieck (MRN, Mathis et al. 1977)-like grain size distribution is adopted for carbonaceous grains and silicates, plus a log-normal distribution for polycyclic aromatic hydrocarbons (PAHs). The photo-electric heating from PAHs is estimated with the prescriptions given in Bakes & Tielens (1994). For every model, the code simulates and outputs detailed thermal (gas and grains) and chemical (with level excitation for selected species) structures of the cloud, as a function of depths, as well as the integrated quantities, such as column densities and line intensities.
For this work, a non-uniform grid of 1367 PDR models has been produced. The explored physical space is the thermal pressure, the intensity of incident UV radiation fields and the depth of the cloud. The thermal pressure, P, ranges from 10 5 to 10 9 K cm −3 . The UV field scaling factor, G UV , ranges from 1 to 10 5 in the Mathis units (Mathis et al. 1983). The visual extinction, A V goes from 1 to 40 magnitudes. Other parameters, common to all models are presented in Table 3 Notes.

Modeling strategy
Alongside the transitions listed in Table 2 Poglitsch et al. 2010) as well as the 13 CO transitions covered by the Herschel SPIRE/FTS bandwidth (J = 5−4 to J = 14−13) are also available. The Meudon PDR code can simulate lineintensities of all these transitions; that is, a priori, all of them could be used to find the best model that matches the observations. Nevertheless, we use only the CO and C lines listed in Table 2 to constrain the physical parameters based on the following three reasons. First, it is not yet clear whether the atomic emissions associated with C + and O share the same origin with the observed CO and C transitions. A comparison of line profiles, resolved both spectrally and spatially, of [CII] 158 µm and CO transitions (J = 2−1, J = 6−5, and J = 13−12) in M17 SW shows that the emissions of C + and CO are hardly associated in their distributions (Pérez-Beaupuits et al. 2012). Spectrally resolved observations of G5.89−0.39, a massive star-forming region in our Galaxy, have also shown a difference between the velocity distributions of [OI] 63 µm and the CO transitions (J = 7−6 and J = 6−5; Leurini et al.2015).
Second, as discussed in Liseau et al. (2006), to match well the model-predicted and observed intensity ratio of [OI] 63 µm/[OI] 145 µm, one has to resort to either extremely high optical depths (A V > 100 mag) or take into account the absorption from the source cloud itself and/or the clouds along the line of sight. With the spectrally resolved [OI] 63 µm line enabled by SOFIA/GREAT (Heyminck et al. 2012), it has indeed been observed in G5.89−0.39 that the emission of [OI] 63 µm is severely contaminated by self-absorption and absorption by clouds along the line of sight (Leurini et al. 2015). Although observational evidence is still scarce at the moment of writing, we judge it best to exclude [CII] 158 µm and [OI] 63 µm when searching for the best physical parameters. As for the [OI] 145 µm emission, since the observation of this transition covers only a small portion (∼20%) of the FOV, we opt to exclude this line as a direct constraint to the models for consistency. However, a good comparison between the model-predicted and observed line intensities of [OI] 145 µm is ensured in our final results, wherever possible.
Third, in order to properly simulate the 13 CO emissions, one must consider the chemical effects due to mutual shielding by its isotopolog, 12 CO. Theoretical computation has indicated that the 13 CO lines can be more efficiently shielded by 12 CO than by itself (Visser et al. 2009). The Meudon PDR code is capable of simulating this process by solving the exact radiative transfer as described in Goicoechea & Le Bourlot (2007). However, this heavy computation would cost days in order to successfully simulate one model, that is, one combination of P, G UV and A V . To efficiently produce the entire model-grid required for the analysis, we adopt the FGK approximation (Federman et al. 1979) that takes only self-shieldings into account. This simplification would result in the underestimation of 13 CO line intensities. Therefore they are excluded from constraining the best physical parameters.
Prior to the search for the best physical parameters on a given pixel based on the nonuniform grid of 1367 PDR models described above, we construct a numerical function that linearly interpolates the grid between its main free parameters, P, G UV and A V , to convert the grid into a continuous 3D grid function. We use radial-basis functions (scipy.interpolate.Rbf) to produce the grid function, M(P, G UV , A V ). In addition to the three free parameters, P, G UV , and A V , an additional free parameter, named scaling factor (φ), is introduced to describe several points. First, this parameter contains the information about the ratio of angular sizes of the observed clouds and beam, that is, φ = Ω cl /Ω beam , where Ω cl is the physical angular size of observed CO clouds and Ω beam is the beam area. Second, it also contains the information about the inclination angle of PDR components in the FOV with the line-of-sight. The Carina I-E region is known to host an edge-on PDR (Brooks et al. 2003), but all detected CO emission in our FOV may not come from edge-on PDRs. So we use face-on line intensities in our grid of PDR models and scale them with φ to consider limb-brightening effect. Third, this parameter also accounts for the possibility that several PDR components are superposed along the line of sight. At each pixel, we search for the combination of φ, P, G UV , and A V , that best matches the observed 12 CO line intensities by minimizing a normalized χ 2 : In Eq. (1), N obs and N p represent the available numbers of the observed CO and C lines and the free parameters (N p = 4), and I i , σ i , and M i (P, G UV , A V ) denote the observed line intensity, its uncertainty, and the simulated line-intensity for the ith transition. As the S/N of the observed CO transitions generally decreases from CO J = 11−10 to CO J = 12−11 and to CO J = 13−12, we rule out the possibility that the intensities for high-J CO transitions beyond the bandwidth of Herschel SPIRE/FTS would be higher than the observed CO J = 13−12 transitions. An artificial penalty in the minimization procedure is forced to better constrain the best-fit solutions. When calculating χ 2 , a penalty is assigned to a given model that predicts I(CO J = 13−12) > I(CO J = 12−11) and I(CO J = 12−11) > I(CO J = 11−10).
We adopt the limited-memory BFGS (L-BFGS) algorithm in our minimization procedure (scipy.optimize.minimize) to search for the local minimum in the χ 2 -space. To ensure that the returned solution of (φ, P, G UV , A V ) indeed corresponds to the global minimum in the χ 2 -space, we randomly choose 100 sets of initial free parameters and output all the solutions found by scipy.optimize.minimize. For a given pixel, the solution that corresponds to minimum χ 2 is chosen to be the final solution. The uncertainties for the constrained solutions are then estimated with 300 Monte-Carlo experiments with an assumption of Gaussian distributions for the uncertainties of the integrated line-intensity for each transition.

CI and CO emission
The CO emission is observed from J = 4−3 to J = 13−12, along with [CI] 370 µm and [CI] 609 µm across the whole FOV (see figures in the Appendix). In the Car I-E region at the northwest, the intensities of the CO and [CI] lines are generally higher by one order of magnitude than those observed in Car I-S and Car I/II. Beyond J = 10−9, the emission is concentrated in the Car I-E and locally in the Car I/II regions. To get a global view on the spatial distribution of the CO excitation over the FOV, we first approximate the excitation temperature (T ex ) of CO transitions observed by the Herschel SPIRE/FTS under the LTE assumption. The reduced-χ 2 for fitting each observed CO SLED with LTE is between 0.5 and 2.5, with an average value of ∼1.5 throughout the FOV. In the Car I-E region, the T ex values are between 70 and 100 K. In the Car I-S region, the T ex values are between 45 and 55 K. In the Car I/II region, the T ex values are generally between 65 and 70 K, except a local (∼0.5 pc in diameter) observation of T ex ∼ 80 K. The observed ratio of [CI] 370 µm/[CI] 609 µmis between 2.9 and 3.7 in the Car I-E region, corresponding to T ex values between 33 and 41 K, which is consistent with the NAN-TEN2 observations (Kramer et al. 2008). In the Car I/II region, the observed ratio of[CI] 370 µm/[CI] 609 µmis around 2.5, corresponding to T ex ∼ 28 K. In the Car I-S region, the [CI] 370 µm and [CI] 609 µm are not detected above a 3σ threshold. The observed spatial trends of the [CI] and CO excitation are generally consistent with each other, showing molecular gas of higher-excitation states in the Car I-E region.
As the sensitivity of the instrument is at the highest around the J = 7−6 transition (Table 2), we present an overview on the morphology of CO emission based on the CO J = 7−6 line intensity. Figure 4 shows the CO J = 7−6 line intensity distribution within our FOV. Comparing with the background optical image shown in Fig. 2, the observed CO gas fills up most of the optically thick region. The most prominent CO and C-emission is found in the Car I-E region, with an average value of 5.3 × 10 −8 W m −2 sr −1 for the CO J = 7−6 intensity, and peaks at the northwest corner of the FOV with a value of 9.3 × 10 −8 W m −2 sr −1 . All 77 pixels (∼2 or ∼1.3 pc in diameter) covered in this region have observed CO 7 − 6 line intensities above the 3σ threshold. This region sits at a projected distance of ∼2.3 pc from Trumpler 14.
Moving to the Car I-S region, the area (∼4 or ∼2.5 pc in diameter) is mostly dominated by the HII region ionized by Trumpler 14. About 30% of the Car I-S region has no observed CO J = 7−6 emission at above the 3σ threshold, while the CO emission observed at above 3σ is distributed along the dust-lane. One double-lined spectroscopic binary (SB2) system, ALS 15 204 (O7.5Vz+O9: V), and a B0 star, ALS 15 203, which are identified as members of Trumpler 14 are found near the center of Car I-S (Apellániz et al. 2016). Along with other stellar members of Trumpler 14, these two stars can have important contribution to the excitation of CO gas in the region. The CO J = 7−6 emission detected in this region is the weakest within our FOV. Among the pixels with detected CO J = 7−6 emission (>3σ), the average intensity is 7.6 × 10 −9 W m −2 sr −1 , which is about an order of magnitude fainter than that in Car I-E.
The Car I/II region at the southeast corner of our FOV is sitting at the intersection of two interstellar clouds, Car I and Car II, as defined by the radio thermal continuum observed at 5 and 0.843 GHz (Gardner & Morimoto 1968;Whiteoak 1994). About 70% of the Car I/II area has CO J = 7−6 emission detected at above the 3σ threshold, and the detected CO J = 7−6 emission has an average intensity of 9.8 × 10 −9 W m −2 sr −1 . This region sits right in the middle between the two most massive star-clusters of the Carina Nebula, Trumpler 14 and Trumpler 16, both of which can provide UV-photons for the observed CO emission.

Physical conditions
We determine the physical conditions (P, G UV , A V and φ) at each pixel following the method described in Sect. 3.2. CO mid-J and high-J lines are good tracers of the gas pressure as it can be seen on Fig. 6 and we expect this parameter, as well as G UV that controls the energy input and photo-destruction rates be well constrained. As only 12 CO lines with J up > 3 and C lines, which  probe only the surface layer of molecular gas are used in constraining the solutions, the total A V parameter of our models is not reliably constrained. So, we do not discuss this parameter in the following. A good agreement between PDR models and observed CO and [CI] line intensities is found at all pixels. To illustrate the general goodness of solutions, comparisons between the observed CO line intensities and those simulated by the Meudon PDR code are presented with the integration of CO lines on the whole FOV (Fig. 7), and with three typical pixels from the Car I-E, Car I-S, and Car I/II regions (Fig. 8). The reduced-χ 2 values vary between 0.2 and 0.7 across the FOV. CO and [CI] line intensities predicted by the Meudon PDR models are generally consistent with the Herschel SPIRE/FTS observations. We also examine whether or not the observed intensities of 13 CO, [OI] 63, [OI] 145 µm, and [CII] 158 µm computed by the PDR models but not used in the minimization procedure are consistent with the observations (see arguments in Sect. 3.2). As our constraints, the 12 CO and [CI] emissions, may not probe the full the Car I-E region, which is dominated by the PDRs, the observed [CII] 158 µm intensities are between 3 and 10 times larger than the model-predicted values. In Car I-S and Car I/II regions, where PDRs no longer dominate the FOV, the observed [CII] 158 µm intensities are about 10 to 50 times larger than the model-predicted values. Generally, this spatial trend is consistent with the results found in Oberst et al. (2011). As the [CII] 158 µm emission can also originate from ionized gas and neutral gas at more diffuse states than the modeled molecular gas based on the CO (Liseau et al. 2006), which leads to the generally over-estimated [OI] 63 µm emission with PDR models. Within our FOV, the ratio of modelpredicted and observed intensities is always larger than 10. The model-predicted [OI] 145 µm transitions, on the other hand, are generally reasonable, and the ratio of model-predicted and observed line intensities are between 1 and 3 over the FOV. Figure 5 shows the maps of derived P and G UV . As each pixel is modeled independently in the process, the smooth appearance on the P and G UV maps reassures that our solutions are not strongly degenerated. Based on the solutions, Car I-E is dominated by high pressure gas (P ∼ 1.5 − 2.2 × 10 8 K cm −3 ) heated by strong UV radiation fields (G UV ∼ 3 × 10 4 ). In the Car I-S region, lower P (≤7 × 10 7 K cm −3 ) and G UV (3000 < G UV < 1.2 × 10 4 ) values are found. The pressure and the intensity of the radiation field then slightly rise into the Car I/II region (P ∼ 10 8 K cm −3 and G UV ∼ 10 4 ). The φ parameter found by our procedure varies between 0.2 and 1.4 in the FOV.
The higher values are found in the Car I-E region. Since in our procedure we use face-on values for simulated line intensities, this is consistent with a PDR seen edge-on that fills the beam as mentioned by previous authors (Brooks et al. 2003;Kramer et al. 2008). The lower values are found in the Car I-S area, consistent with a region where molecular gas partially fills up the beam.
The G UV values obtained by our comparison between PDR models and observations are larger than the ones cited in previous works. In particular, Brooks et al. (2003) estimated the FUV flux in several ways. First, based on the Kaufman et al. (1999) PDR models, they find a FUV scaling factor between 600 and 10 4 in the Habing unit. They also considered the most massive stars in the Tr 14 cluster and estimated the resulting flux at the [CII] peak location (10:43:22 -59:34:45), finding 1.4 × 10 4 in the Habing unit, that is, 10 4 in the Mathis units.
We re-investigate the FUV luminosity that illuminates the molecular gas observed in the FOV. The spatial variation of our derived G UV across the FOV suggests that the CO clouds may be distributed at different distances from the illuminating stars. As the massive stellar members of Trumpler 14 and 16 are well identified (Vazquez et al. 1996;Sana et al. 2010;Hur et al. 2012;Sota et al. 2014;Apellániz et al. 2016;Alexander et al. 2016), we are able to compare the derived G UV with the estimated UVphotons input by the nearby OB-stars. Table 5 lists all observed OB-stars within a projected distance of 7 pc (∼7.5 ) from the center of our FOV with updated FUV luminosity. If one is unaware of the spatial distributions of stars and naively assumes that all the listed stars are located at the center of their host star-clusters, and that the observed CO gas lies at the same plane perpendicular to our line of sight at a 2.3 kpc distance, the estimated stellar G UV values and the G UV values derived from the Meudon PDR models agree poorly with each other (left panel in Fig. 9). However, when the sky coordinates of the nearby OB-stars are taken into account (still coplanar with the CO gas), the comparison between the estimated stellar G UV values and the G UV values derived from the Meudon PDR models slightly improves in the Car I-E and Car I/II regions (center panel in Fig. 9).
Two stellar members of Trumpler 14, an O-star binary system, ALS 15 204, and a B0 star, ALS 15 203, are observed in the Car I-S region (near the pointing #3 in Fig. 2). Observations with the 4 m Anglo-Australian Telescope (AAT) have suggested that a separate subcluster may be formed around these two stars (Alexander et al. 2016). Judging from the spectrally resolved CO J = 1−0 observations 8 (Rebolledo et al. 2016) across the six central pointings (see Fig. 2) from our FOV in Fig. 10, the CO 1 − 0 emission at Car I-S region is dominated Left panel: all OB-stars are assumed to be at the center of their host star-clusters, i.e., Trumpler 14 and 16. The observed CO gas is assumed to lie on the plane perpendicular to the line of sight at a 2.35 kpc distance. Center panel: distribution of OB-stars is considered coplanar at their observed (α, δ) coordinates. The geometry of the CO gas is assumed to be coplanar as in the textitleft panel. Right panel: distribution of OB-stars is considered coplanar at their observed (α, δ) coordinates. The geometry of the CO gas is assumed as the shape indicated in Fig. 11. by the velocity component at around −13 km s −1 (pointings # 3, # 4, and # 5). Assuming that the subcluster carves out a cavity of CO gas in the Car I-S region, the observed CO emission from the region may be associated to the cloud sitting further away from the observer along the line of sight. With a simple cloud geometry demonstrated in Fig. 11, which is created by moving the cloud directly over ALS 15 204 and ALS 15 203 in Car I-S further by ∼3 pc 9 , ensuring the spatial continuity of the cloud along constant slopes into Car I-E and Car I/II, and taking the coordinates of nearby OB-stars into account, the G UV values derived with the Meudon PDR code come to an excellent agreement with the expected G UV contributed by the OB-stars in the field without scaling (right panel in Fig. 9). Even if we assume that the distance to Trumpler 14 is ∼2.9 kpc (Hur et al. 2012), in which case the expected G UV contributed by the OB-stars would decrease by a factor of ∼1.5, the agreement still holds within the estimated uncertainties.
In our simple and crude geometry assumption, all the stars in the vicinity are located on the same plane, which is unlikely the reality. Unfortunately, given the uncertainties in our observations and models, we do not have sufficient information to accurately constrain the positions of stars and clouds along our line of sight. However, the excellent agreement displayed in Fig. 9 demonstrates that our simple assumption may be a good approximation of reality.

UV radiation field in Carina
The CO emitting gas observed in this work only covers a small area (∼2 × 5 pc 2 ) of the Carina Nebula (Fig. 1). In this subsection, we discuss our findings under the frame of the whole Carina Nebula. The Carina Nebula is famous for its bipolar structure with its major axis almost perpendicular to the Galactic plane. In its northern part, expansion along the Galactic plane caused by the formation of young star clusters is inhibited by the dense molecular gas (Smith et al. 2000). One of the biggest mysteries in the Carina Nebula is its supernova history.
Based on the observed line-splitting of ionized gas and the total gas mass estimated from the cool dust mass, A&A 618, A53 (2018)   However, the luminosity for this star is uncertain and can possibly be five times brighter (Hamann et al. 2006). Smith & Brooks (2007) estimated the total kinetic energy of the nebula to be 8 × 10 51 erg, which represents 30% of the total mechanical energy input by the stellar wind during the past 3 Myr (Smith 2006a). The lack of radio or X-ray synchrotron emission indicates that a recent supernova activity is unlikely. However, confirmation of the diffuse X-ray emission which offsets from the young massive star clusters and the discovery of an isolated neutron star, 2XMM J104608.7−594306, have led to a possible supernova history in the Carina Nebula (Townsley et al. 2011;Pires et al. 2015). In the vicinity of our FOV, an observed B-type star (Tr 14-29 in Table 5, also known as Trumpler 14-MJ 218) is moving with a high velocity (∼100 km s −1 ) toward the center of Trumpler 14, directly away from 2XMM J104608.7−594306. If Tr 14-29 is the runaway star associated with 2XMM J104608.7−594306 during a past supernova explosion, its estimated flight time is ∼(1.1 − 3) × 10 4 yr (Pires et al. 2015). In an alternative scenario where Tr 14-29 originates from the open star cluster, Trumpler 16, its flight time would be 10 5 yr (Ngoumou et al. 2013). Both estimated flight times are roughly consistent with the absence of a supernova remnant in the Carina Nebula. All these characteristics suggest that the region in our FOV may be associated with a Galactic superbubble. The observed physical conditions of PDRs can thus provide a suitable analog for the observations of more faraway starburst regions hosted in similar interstellar environments, such as the famous 30 Doradus.
Comparing the CO SLEDs from the Car I-E, Car I-S, and Car I/II regions (see Fig. 8), it is clear that the Car I-E region is Fig. 11. Estimated 3D cloud structure in this work (gray mesh). The blue star symbols indicate the OB-stars in the FOV with their sizes in proportion to the FUV luminosities of stars. The Z-axis indicates the offset along the line of sight where Z = 0.0 is at 2.3 kpc from us. Stars that are outside the FOV (Fig. 2) toward the south are displayed in fainter colors. dominated by CO gas of high thermal pressure and strong radiation fields (P = 1.8 × 10 8 K cm −3 and G UV = 1.7 × 10 4 ). The thermal pressure and the strength of radiation decrease (P = 5.6 × 10 7 K cm −3 and G UV = 3.4 × 10 3 ) in the Car I-S region and then increase again in the Car I/II region (P = 1.1 × 10 8 K cm −3 and G UV = 1.5 × 10 4 ). Within our FOV, the PDRs in the Car I-E region have been previously studied by several groups. Based on the approximate stellar composition and the FIR continuum observed in Car I, Brooks et al. (2003) and Mizutani et al. (2004) (Oberst et al. 2011). With a clumpy PDR model, using the observed CO J = 2−1, J = 4−3, J = 7−6, and the two atomic carbon lines, [CI] 370 and 609 µm, the strength of the radiation field in Car I-E is estimated to be ∼3200 in the Habing unit (∼2500 in the Mathis unit) (Kramer et al. 2008). In a big picture, these results indicate that the UV-photons and stellar winds provided by nearby massive stars are more than sufficient in supporting the global heating of CO and [CI] emissions observed in this region.
Using the Meudon PDR models, the G UV values constrained from the Herschel SPIRE/FTS observations in the Car I-E region range between 10 4 and 3 × 10 4 , which are generally higher than the values found in the literature. However, the excellent agreement of the model-derived G UV values with the estimated radiation fields from the nearby stellar composition appears throughout our FOV, as demonstrated in Fig. 9. Our analysis has confirmed that the chief energy source for the observed CO emissions is the UV-photons input by the massive stars in the vicinity. We would like to note from our results that in the Car I-E region, the model-derived G UV values tend to be marginally higher than the G UV contributed solely by the OB-stars. The emission of CO observed along the ionization front may also be supported by local heating sources. In Car I-E, narrow-band imaging-observations from the optical to the mid-infrared have revealed that many massive (B4 to A0) young stellar objects (YSOs) are distributed on or just behind the ionization front (Tapia et al. 2006). Along the ionization front, the modelderived thermal pressure values are the highest in our FOV (∼2 × 10 8 km cm −3 , see Fig. 5). Interestingly, along the high-pressure ionization front, at the very edge of our FOV in the northwest, the highest values for the CO and [CI] emissions are observed (see figures in the Appendix). At this location (10:43:23.25, -59:33:56.9), a massive class I YSO, V723 Car, has been observed to go into outburst just before 2003 (Tapia et al. 2015). All 12 CO transitions covered by the Herschel SPIRE/FTS and 13 CO, up to J = 8−7, is unambiguously detected here. This result suggests a potential contribution of local heating sources. Unfortunately, we cannot derive the physical conditions from this location due to the variation of instrument resolution, but it will be an interesting target for future high-resolution observation for understanding the interplay between the formation of massive YSOs and their host environment.
The location of Car I/II region has almost the same distances to the centers of Trumpler 14 and 16 and lies at the intersection of Car I and Car II. Based on the projected distances from the O-and B-stars in the vicinity (see Table 5), the G UV -values at the Car I/II region should be comparable to or weaker than those at the Car I-E region. This is the case for most of the pixels in the region. Globally, in Car I/II, the dust mass and temperature derived from the photometry data in our FOV are quite low, however CO is unambiguously observed (Wu et al., in prep.). This result is in agreement with the comparison of dust-mass estimated with the continuum at 1.2 mm (9 M , Brooks et al. 2005), and gas-mass estimated with CO J = 2−1 (91 M (LTE) or 411 M (virial), Rathborne et al. 2002). The large discrepancy between the gas-mass estimated under LTE and virial assumption implies that the observed cloud is not gravitationally bound but supported by external radiation fields (Rathborne et al. 2002).
One thing that has been puzzling from our results is the local high G UV values (∼5 × 10 4 ) derived in the Car I/II region of our FOV (see Fig. 5). At this location (10:44:03.89, -59:37:57.63), CO emission is well detected up to J = 12−11. Under the LTE assumption, the T ex value derived locally is ∼80 K, while T ex is generally below 70 K in the Car I/II region. This point also appears as the only outlier in Fig. 9, hinting that the UV-photons from nearby OB-stars may not be able to support the CO excitation. Observations in CO J = 1−0 with Mopra show a velocity component centered at −5 kms −1 , which appears only at this point in our FOV (see dashed-line in Fig. 10), implying that the structure of CO gas may differ locally. The CO cloud observed here is likely heated partially by stellar members in Trumpler 14 and partially by the two nearby early-type stars, HD 93162 and ALS 15 210, of Trumpler 16. ALS 15 210 has been classified as an O 3.5 I star. The Wolf-Rayet star, HD 93162, also known as WR 25, is possibly the most luminous known star in our Galaxy (up to L = 6.3 × 10 6 L , Hamann et al. 2006), but it is still not sufficient to account for the local high G UV value. The UV-photons input by the nearby massive stars simply cannot explain the local model-derived G UV values. Based on the modeling results of the 1−24 µm spectral energy distribution in the Carina Nebula, a type 0/I YSO of intermediate-mass with hard X-ray signature has been detected at (10:44:03.086, -59:37:48.41) (Povich et al. 2011) and listed in the Pan-Carina Young Stellar Object Catalog (PCYC) as PCYC 419. However, at the physical scale of ∼0.2 pc, we do not expect the outflow from an intermediate-mass YSO to dominate the heating of CO gas. To further understand the local heating mechanisms for the observed CO gas, high-resolution observations from the region are necessary.

CO at PDR interfaces
For more than 30 yr, PDRs have been studied to understand the transition between atomic and molecular gas as well as the impact of the radiative feedback of massive stars on their parent cloud (Tielens & Hollenbach 1985;Sternberg & Dalgarno 1989). Before Herschel observations, the most accessible lines probing the very edge of PDRs and the physics and chemistry that take place in these interfaces were the atomic lines of O, C + and C as well as molecular lines of H 2 and a few CO lines. Comparisons of these observations to PDR models lead to the conclusion that the edge of PDRs can be seen as a clumpy medium in which dense clumps are embedded in a more diffuse medium (Stutzki & Guesten 1990;Parmar et al. 1991;Meixner & Tielens 1993;Tauber et al. 1994;Hogerheijde et al. 1995;Andree-Labsch et al. 2017). In this scenario, UV photons can penetrate deeply into the cloud through the interclump medium (n H ∼ 10 4 −10 5 cm −3 ) and heat the edge of the dense clumps (n H ∼ 10 6 −10 7 cm −3 ) where H 2 vibrational lines and CO mid and high-J lines are produced. More recently, Herschel observations of CO rotational emission in protostars, PDRs and extragalactic regions also showed high-J excitation. Shapes of CO SLEDs in different environments have been compared by Indriolo et al. (2017). In protostellar systems such as Orion KL, the very high-J excitation (J up ∼ 40) indicates energy inputs through outflows and shocks (Goicoechea et al. 2015). In extragalactic regions, the sources of energy responsible for CO excitation are less clear because of the mixture of several emitting sources in the beam. Rosenberg et al. (2015) suggest that CO excitation in some Class I galaxies may be dominated by UV heating whereas Class III galaxies would require mechanical heating (shocks and turbulence) or X-rays heating.
Here, we find that highly illuminated (G UV ∼ 10 4 ) PDR models with high thermal pressures (P ∼ 10 8 K cm −3 ) can successfully explain CO excitation in the Carina Nebula without introducing a clumpy cloud structure. UV photons, via the photo-electric effect on grains, are the main heating mechanism. Isobaric models have been proposed as an alternative to the classical clump approach for the Orion Bar (Marconi et al. 1998;Allers et al. 2005). An isobaric model is found to successfully explain the CH + emission in the Orion Bar (Nagy et al. 2013). More recently, it is suggested that the PDRs can support several tens of atomic and molecular lines including H 2 lines and CO lines from J = 4−3 to J = 23−22 in two prototypical PDRs, the Orion Bar and NGC 7023 NW (Joblin et al. 2018).
It is suggested that hot chemistry takes place at the edge of PDRs and is responsible for the mid-J and high-J CO lines that are observed (Joblin et al. 2018). In our models, we simulate the formation of H 2 with the formalism described in Le Bourlot et al. (2012) that considers the Langmuir-Hinshelwood and Eley-Rideal mechanisms. This detailed treatment increases the formation rate of H 2 compared to the classical value determined in diffuse interstellar gas by Copernicus and FUSE, 3 × 10 17 cm 3 s −1 (Jura 1974;Gry et al. 2002) that is often used in astrochemical models. Typically, in the Car I-E area, we find a three to four times higher formation rate at the edge of the PDR. This enhanced formation rate allows H 2 to form close to the edge of the PDR in warm gas (a few hundred Kelvins depending on physical conditions). In this warm gas, some key chemical reactions with H 2 , usually inhibited because of large activation energies, can be efficient because of first, the relative velocity between the reactants and second, the internal energy stocked in H 2 rovibrational levels. That is the case for the C + + H 2 reaction that forms CH + with an activation energy of 4537 K. We compute the rate of this reaction with the formalism of Agúndez et al. (2010).
Several chemical routes that start with CH + lead to the formation of CO. As a consequence, in bright PDRs, the formation of CO starts as soon as H 2 self-shielding comes into effect. Figure 12 presents the density profiles of H, H 2 , C + , C and CO (total and in its levels J = 4, 8 and 13) as a function of A V for a model with P = 3.5 × 10 8 K cm −3 and G UV = 2 × 10 4 . CO level J = 13 is populated before the C/CO transition. Therefore, a fraction of the emission of CO in mid-J and high-J levels comes from a slab of the PDR in front of the C/CO transition. Recent ALMA observations of the Orion Bar presented in Goicoechea et al. (2016Goicoechea et al. ( , 2017 show the same trend of chemical stratification. In particular, SH + emission whose formation path is similar to the one of CH + but with an activation energy threshold of 9860 K, shows a clear front close to the edge of Bar.  Figure 13 presents the thermal pressure as a function of G UV for all the pixels of our Carina FOV. As the pressure, P, is derived from CO and C lines in this work, it is generally associated with neutral gas. Except in the Car I-S region, which is not very bright in the CO emission and dominated by the HII region, we fit a power-law relationship between P and G UV values in Fig. 13. In the Car I/II region, some pixels seem to have the same thermal pressure for different G UV while for other ones P increases with G UV . In the Car I-E region, which is dominated by the bright PDR, the relationship is the most prominent. From our observations, the empirical relationship between P and G UV is:

P-G UV relation
In the extragalactic observation, a similar relationship has been suggested in Wu et al. (2015). In the Galactic PDRs, a similar relationship between P and G UV values spanning over more than three orders of magnitude has also been observed (Joblin et al. 2018). The observed relationship suggests that the UV radiation field of massive stars may create a high pressure layer at the edge of PDRs and the physical mechanism at play could be the UV induced photo-evaporation as described in Bertoldi (1989), Bertoldi & Draine (1996). As the Meudon PDR model describes a stationary system, a hydrodynamical PDR code is required to understand the exact mechanism leading to the compression of the gas and to quantify its effect ).

Conclusion
The Herschel SPIRE FTS uses a pioneering design to simultaneously observe a broad range of CO transitions in spectral imaging mode which enables the derivation of physical properties of molecular gas in a macroscopic field of view. Combined with the state-of-the-art Meudon PDR models, we obtain the following results: -The CO emission is observed in the Carina Nebula, near the young star cluster, Trumpler 14, from J = 4−3 to J = 13−12 across the Car I-E, Car I-S and Car I/II regions ( Fig. 2 and figures in the Appendix). In the Car I-E region at the northwest of our FOV, the intensities of the observed CO are generally higher by one order of magnitude than those observed in Car I-S and Car I/II. Beyond J = 10−9, the emission is concentrated in the Car I-E and locally in the Car I/II regions, hinting at higher excitation temperatures of CO gas in these regions. -We successfully account for the observed CO and [CI] emission with high-pressure and high-FUV-illumination PDR models. Comparing with previous PDR modeling results found in the same region, the G UV values derived in this work are generally higher (Brooks et al. 2003;Mizutani et al. 2004;Kramer et al. 2008;Oberst et al. 2011). However, the excellent relationship observed between the model-derived G UV values and the estimated radiation fields (Fig. 9)  observations. -The thermal pressure derived in this work is of the order of magnitude 10 8 K cm −3 across our FOV, with the highest values found at the dissociation front in the Car I-E region. Although the spatial resolution of our observations prevents us from obtaining a detailed cloud structure at the dissociation front, our results are in general agreement with the recent ALMA observations of the Orion Bar that the observed cloud edge is compressed by a high-pressure (2 × 10 8 K cm −3 ) wave propagating into the molecular cloud (Goicoechea et al. 2016. -We find a relation between the thermal pressure at the edge of the PDR and the intensity of the FUV radiation field, P ∼ 2 × 10 4 G UV . This relation is similar to the one found by Wu et al. (2015) and Joblin et al. (2018). It suggests that massive stars may create a high-pressure slab at the edge of PDRs. -Two local (∼0.2 × 0.2 pc) regions in our FOV deserve further investigation with higher-resolution observations: (1) At the very edge of our FOV (10:43:23.25, -59:33:56.9), a massive class I YSO, V723 Car has been observed to go into outburst just before 2003 (Tapia et al. 2015). All 12 CO transitions covered by the Herschel SPIRE/FTS and 13 CO, up to J = 8−7, are unambiguously detected here. Unfortunately, we cannot derive the physical conditions from this location due to the variation of instrument resolution, but it will be an interesting target for understanding the interplay between the formation of massive YSOs and their host environment.
(2) In the Car I/II region, we find high G UV values (∼5 × 10 4 ) locally at (10:44:03.89, -59:37:57.63). At this position, CO emission is well detected up to J = 12−11, and the LTE-approximated T ex value also implies locally excited gas. We cannot identify the possible heating source(s) based on the available observations. To further understand the local heating mechanisms at these locations, high-resolution observations from the region are necessary.

SPIRE FTS in M83
We show the observed maps of the lines listed in Table 2 by the SPIRE FTS. All the maps are displayed in the observed spatial