Constraining the geometry of the reflection nebula NGC 2023 with [O i]: emission & absorption

ABSTRACT

We have mapped the NGC 2023 reflection nebula in the 63 and 145 |$\mu$|m transitions of [O i] and the 158 |$\mu$|m [C ii] spectral lines using the heterodyne receiver upGREAT on SOFIA. The observations were used to identify the diffuse and dense components of the photon-dominated region (PDR) traced by the [C ii] and [O i] emission, respectively. The velocity-resolved observations reveal the presence of a significant column of low-excitation atomic oxygen, seen in absorption in the [O i] 63 |$\mu$|m spectra, amounting to about 20–60 per cent of the oxygen column seen in emission in the [O i] 145 |$\mu$|m spectra. Some self-absorption is also seen in [C ii], but for the most part it is hardly noticeable. The [C ii] and [O i] 63 |$\mu$|m spectra show strong red- and blue-shifted wings due to photoevaporation flows especially in the south-eastern and southern part of the reflection nebula, where comparison with the mid- and high-J CO emission indicates that the C⁺ region is expanding into a dense molecular cloud. Using a two-slab toy model the large-scale self-absorption seen in [O i] 63 |$\mu$|m is readily explained as originating in foreground low-excitation gas associated with the source. Similar columns have also been observed recently in other Galactic PDRs. These results have two implications: for the velocity-unresolved extragalactic observations this could impact the use of [O i] 63 |$\mu$|m as a tracer of massive star formation and secondly, the widespread self-absorption in [O i] 63 |$\mu$|m leads to underestimate of the column density of atomic oxygen derived from this tracer and necessitates the use of alternative indirect methods.

1 INTRODUCTION

The fine-structure lines of [C ii] at 158 |$\mu$|m and [O i] at 63 and 145 |$\mu$|m are the most important cooling lines in the far-infrared (FIR) and have long been used to study photon-dominated regions (PDRs). In the rest of the paper, we refer to the [O i] 63 |$\mu$|m and [O i] 145 |$\mu$|m fine-structure lines simply as [O i] 63 and [O i] 145. PDRs are regions where far-ultraviolet (FUV; 6 eV < hν 13.6 eV) radiation from young massive stars dominates the physics and the chemistry of the interstellar medium and correspond to the transition from ionized to molecular gas (Hollenbach & Tielens 1997). PDRS were first studied extensively with the Infrared Space Observatory, though the lack of spectral and spatial resolution was a severe limitation. This changed with the Kuiper Airborne Observatory (in limited regions) and with the Heterodyne Instrument for the FIR on Herschel, which had both velocity resolution and sensitivity to even enable determination of the optical depth of [C ii] by also observing ¹³[C ii] (Ossenkopf et al. 2013). After Herschel, this work was taken over by first GREAT on SOFIA and later with the array receiver upGREAT. These studies have found that the [C ii] emission often can be moderately optically thick and sometimes significantly self-absorbed (Mookerjea et al. 2018; Graf et al. 2012; Guevara et al. 2020), which may underestimate [C ii] column densities by a factor of two to three. The [O i] 63 |$\mu$|m can also be optically thick and self-absorbed, sometimes even more strongly than [C ii] (see e.g. Mookerjea et al. 2021).

NGC 2023, illuminated by the B2 star HD 37903, is one of the best-studied reflection nebulae. It is nearby, ∼400 pc, and exhibits a strong, nearly edge-on PDR in which the cavity expands into the surrounding dense molecular cloud. It has served as a test bed for exploring PDR models (e.g. Kaufman, Wolfire & Hollenbach 2006). It was mapped with GREAT on SOFIA in [C ii] and CO(11–10) with the single pixel L1/L2 mixer (Sandell et al. 2015), who also mapped the nebula in ¹³CO(3–2), and various transitions of CO from J=3–2 up to J=7–6.

Here, we have additionally mapped the NGC 2023 PDR in the 63 and 145 |$\mu$|m fine-structure transitions of [O i] and at the same time obtained a deeper map of [C ii] using upGREAT on SOFIA. The [O i] observations enable us for the first time to study the physical conditions in the south-eastern and southern part of the nebula, where the C⁺ region is expanding into the dense surrounding molecular cloud as well as estimate the relative contribution of the dense and diffuse PDR gas to the [C ii] emission.

The newly obtained observations are analysed in combination with the previously published maps of (i) ¹³CO(3–2), CO(6–5), CO(7–6) observed with APEX with beamsizes of 18|${_{.}^{\prime\prime}}$|5, 9|${_{.}^{\prime\prime}}$|1, and 7|${_{.}^{\prime\prime}}$|7, respectively and (ii) CO(11–10) observed with SOFIA/GREAT with a beamsize of 23 arcsec (Sandell et al. 2015).

2 OBSERVATIONS

The reflection nebula NGC 2023 was observed on a 90 min leg with the SOFIA upGREAT¹ receiver (Risacher et al. 2016) in GREAT consortium time on SOFIA flight #668 (Project ID: 83_0731) out of Palmdale, CA, on 2020 March 6. Although NGC 2023 was observed at low flight altitude, 11.6 km (38 000 feet), the observing conditions were very good with a precipitable water vapour around 5–10 |$\mu$|m. The upGREAT was in the Low-Frequency Array/High-Frequency Array (LFA/HFA) configuration with both arrays operating in parallel. The V polarization of the LFA was tuned to [C ii] at 1.900 5369 THz while the H polarization was tuned to the [O i] 145 |$\mu$|m line at 2.060 068 86 THz. The HFA was tuned to the [O i] 63 |$\mu$|m at 4744.777 49 THz. We made total power on-the fly (OTF) maps in classical position switched mode. The reference position was at +390″,−70″ relative to HD 37903 (α_J2000: 05^h41^m38.39^s, δ_J2000: −02^○15′32|${_{.}^{\prime\prime}}$|5). The OTF map was done in two coverages, scanning in both x and y. The map was rotated by 45^○. The sampling was done every 3 arcsec with a sampling time of 0.3 s per dump. This resulted in maps for the LFA array of ∼4|${_{.}^{\prime}}$|9 × 3|${_{.}^{\prime}}$|9, while the map size for the HFA array was 3.0′× 2|${_{.}^{\prime}}$|1, which was enough to cover the SE quadrant of the reflection nebula, where [O i] was expected to be strong.

The observations were reduced and calibrated by the GREAT team. The GREAT team also provided beamsizes (14|${_{.}^{\prime\prime}}$|1 for [C ii], 13|${_{.}^{\prime\prime}}$|0 for [O i] 145 |$\mu$|m, and 6|${_{.}^{\prime\prime}}$|3 for [O i] 63 |$\mu$|m) and beam efficiencies derived from planet observations. The data were corrected for atmospheric extinction, flat fielded and calibrated in T_mb. Further processing of the data was made with the class² software.

We also retrieved some high-quality [C ii] data from the SOFIA archive, which extended the [C ii] map to the north. This data set was from project 02_0090 (PI. Els Peeters) and fully calibrated in the SOFIA archive. All data in this paper are presented in main-beam temperature scale.

3 RESULTS

3.1 Morphology & kinematics of the region

Sandell et al. (2015) used maps of multiple CO rotational transitions combined with the [C ii] fine-structure line at 158 |$\mu$|m to probe the overall morphology of the reflection nebula NGC 2023. These observations showed that the [C ii] emission traces an expanding ellipsoidal shell-like PDR region region, which is surrounded by a hot molecular shell. In the south-east, where the [C ii] region expands into a dense, clumpy molecular cloud ridge, they also detected emission from high-J CO lines, apparently originating in a thin, hot molecular shell surrounding the [C ii] emission. These authors found that there was a clear velocity gradient across the nebula, with the emission being more blue-shifted in the south and south-east and more red-shifted in the north and north-west, some of which appeared to be due to an expansion of the nebula. Sandell et al. (2015) also noted that high-angular resolution images of vibrationally excited H₂ and PAH emission suggested that the PDR was far from smooth, exhibiting lumpy ridges, and bright filamentary structures. The [C ii] and high-J CO maps looked smoother, partly because they had insufficient spatial resolution to see such details. Here we have additionally mapped the fine-structure transitions of [O i] at 63 and 145 |$\mu$|m, which are good tracers of the warm and dense regions of the PDR. The upper energy levels of the 63 |$\mu$|m and 145 |$\mu$|m [O i] transitions are 227.7 K and 326.6 K respectively, while formally the critical densities for these transitions are 5 10⁵ cm⁻³ and 5 10⁶ cm⁻³ respectively, for collisions with H₂ (Goldsmith 2019). The [O i] 63 emission, however, can be very optically thick and self-absorbed, which was seen toward the S1 PDR in RhoOph (Mookerjea et al. 2021) and the same is true for NGC 2023. However, here we have also observed [O i] 145, which is unaffected by self-absorption. The [O i] 145 emission is largely optically thin and shows that the dense PDR has a lumpy and filamentary structure. The [O i] 145 map looks qualitatively similar to the [C ii] map (Fig. 1), but the [O i] emission is concentrated to the dense south-eastern and southern part of the PDR and more structured, while the [C ii] emission is rather smooth. In the north and north-west the gas densities are too low to collisionally excite [O i]. The [O i] 63 map, which has a spatial resolution of 6|${_{.}^{\prime\prime}}$|3, shows that the emission is extremely patchy over the area where it is detected (Fig. 2). The emission is essentially completely self-absorbed in the south-east and south, while it is seen in emission in the east and north-east, where the absorption from the foreground cloud is less severe. Even here it is strongly affected by absorption from colder foreground gas giving the integrated line emission a very patchy and lumpy appearance.

We also compare the [O i] 145 map with maps of CO(7–6) and CO(11–10) (Sandell et al. 2015), which also trace the hot, dense PDR as well as with C91α, which was mapped by Wyrowski et al. (2000; Fig. 3).

The overall morphology of [O i] 145 matches CO(7–6) and CO(11–10) quite well. Both appear to trace hot gas in the leading edge of the expanding PDR shell. That CO(7–6) emission is not seen from the entire PDR where [O i] is detected, is likely due to the fact that the emission is strong enough only from a narrow shell and needs to be viewed almost tangentially to have enough column density to be detected. The emission from the remaining hot gas is likely below the sensitivity of our observations. C91α shows a different morphology than CO(11–10), although they both appear to originate in approximately the same region. However, where C91α is strong, CO(11–10) is faint, like in the southern part of the PDR shell, where C91α shows two strong peaks. Both of these peaks are seen in [O i] 145, which appear to trace the same gas as C91α. In the east, where C91α is faint there is no clear peak in [O i] 145 either. We note that although the [O i] 145 traces the high-density ridge and the clumps in the ridge, the emission is slightly offset compared to CO(11–10) and C91α, which appears to be further outside in the PDR shell, i.e. closer to the cold molecular cloud (Fig. 3). Based on the detailed analysis of the high-J CO lines, Sandell et al. (2015) estimated the hot molecular shell to be between 90–120 K and have densities n ∼ 10⁵–10⁶ cm⁻³.

The velocity-channel map of [O i] 145 emission clearly show that the dense PDR clumps along the ridge are also at different velocities with the south-western clumps being more red-shifted than the north-eastern clump (Fig. 4). At v = 9, km s⁻¹ the region north-east of HD 37903 dominates the [O i] 145 and [C ii] emission. This emission peak, at ∼+27″, +18″, coincides with the strongest [C ii] emission peak, ∼55 K, in the whole map. It is also seen as a distinct peak in CO(6–5), CO(7–6), and CO(11–10), but it is relatively faint compared to the emission at ∼v = 10.5 km s⁻¹, which dominates the emission in the south-eastern part of the PDR. This suggests that the north-eastern emission peak may have somewhat lower gas densities than the PDR emission in the ridge.

Optical images show that the NGC 2023 reflection nebula is quite large with a size of ∼10 × 10 arcmin. In the near- and mid-IR it is only about half the size, ∼5|${_{.}^{\prime}}$|2 × 5|${_{.}^{\prime}}$|3 (Mookerjea, Sandell & Jarrett 2009). This is the regime where the emission is dominated by the PDR. In the near-IR one can see that HD 37903 illuminates a bright ridge ∼120–150 arcsec to the north, which intercepts some of the FUV radiation. This ridge is dominated by fluorescent H₂ emission, see e.g, the narrowband H₂ imaging by Field et al. (1998). The feature Field et al. (1998) calls the IR triangle is quite bright in [C ii] (Fig. 1). The peak of the [C ii] emission (+27″, +18″) coincides with the head of what Field et al. (1998) refer to as the Seahorse. It is clear that the PDR shell opens up to the north, as well as in a region east of the IR triangle, Fig. 1. Even in the south-east, where the C⁺ region is surrounded by the dense molecular cloud, one can see that FUV radiation escapes, like the tongue of [C ii] emission sticking out at a position angle (measured counterclockwise from North) of ∼110^○ (Fig. 1).

The position–velocity diagrams, Fig. 5, provides us with an alternative way to examine the morphology of the PDR and what the different PDR tracers show us. There is a clear velocity gradient from south to north with the emission in the northern part of the PDR being more red-shifted, seen in all the tracers, even [O i] 145, although the [O i] 145 map only partially covers the PDR ridge. This suggest that the IR triangle and the northern ridge is on the backside of the PDR. Near HD 37903, the [C ii] emission shows both blue- and red-shifted emission, suggesting that the whole C⁺ region is embedded in the cloud, because one sees photoevaporation flows both from the front and the backside of the PDR. The CO(4–3) and CO(6–5) show an outflow just south of HD 37903 and strong blue- and red-shifted emission, where the position–velocity diagram crosses the outflows from Sellgren C and D. Although [C ii] also shows strong blue-shifted emission in the south, it is almost certainly due to strong photoevaporation, as one can see from the velocity-channel maps (Fig. 4). These maps show strong blue-shifted emission over the whole PDR in the south-east. The east–west position–velocity diagram shows the strong blue-shifted photoevaporation flow inside the PDR shell and it appears that some FUV radiation ‘leaks’ through the shell, i.e. one can see [C ii] emission at distances >100 arcsec from HD 37903 at roughly the same velocity, ∼9.9 km s⁻¹. Closer to the star the [C ii] emission appears to be somewhat self-absorbed. [O i] 145 shows blue-shifted emission east of HD 37903, indicating that the photoevaporating gas is so dense that it is seen even in [O i] 145. To the west, the [O i] 145 becomes more red-shifted, which is not really noticeable in [C ii]. The CO lines pick up the blue-shifted outflow slightly west of HD 37903. The position–velocity diagram at PA 110^○ shows one of the regions where [C ii] breaks through the PDR shell. However, it looks about the same as in the east–west cut. Towards the west-north-west the [C ii] splits up into two velocity components. However, both [O i] 145 and ¹³CO(3–2) only show a single velocity component roughly halfway between the two [C ii] peaks, indicating that [C ii] must be self-absorbed.

The position–velocity diagram along the ridge cuts through the two strong C91α peaks. The southern one is faint in [C ii] while it is quite prominent in [O i] 145. Here the emission is shifted toward ∼10.7 km s⁻¹, whereas the clump to the north-east of it is ∼10.2 km s⁻¹ (Fig. 5). The velocity of the [O i] 145 emission appears to follow that of C91α (Wyrowski et al. 2000), while these emission clumps are unnoticeable in [C ii]. The CO lines show strong blue-shifted emission from a molecular outflow powered by the young star Sellgren D.

3.2 Spectral profiles of [O i] lines

Based on the model proposed by (Sandell et al. 2015) the spectral profiles of CO and [C ii] detected NGC 2023 consist of several velocity components: the quiescent cloud seen in low-J CO lines and ¹³CO, the PDR, and line wings from high-velocity molecular outflows seen in the outskirt of the reflection nebula. The [C ii] spectra also show prominent blue- or red-shifted wings from photoevaporation flows. In order to see which of these components contribute to the [O i] emissions, we have selected six positions, some of which already were analysed by (Sandell et al. 2015). All spectra are smoothed to a common angular resolution of 15 arcsec (Fig. 6). We have also fitted the [O i] 145 spectra using Gaussian profiles consisting of one or two velocity components as appropriate and compared the fitted parameters with those derived for [C ii] (this work) and CO(6–5), (7–6), and (11–10) transitions (Table 1). We have not used [O i] 63 |$\mu$|m, because the line is strongly self-absorbed over the whole area that we have mapped. The [O i] 145 spectra are all single-peaked, although we do see clear blue-shifted line wings in several positions, like at 40″, −40″ and 23″, −60″, where [C ii] shows very strong blue-shifted wings. The position 23″, −60″, which coincides with the C91α emission peak 2, also shows blue-shifted emission in [O i] 145, and even in C91α. At position 60″, −60″, the [C ii] line is strongly self-absorbed at the velocity of [O i] 145, making the blue-shifted line wing appear stronger than the emission from the PDR (Table 1). Sandell et al. (2015) interpreted the double-peaked [C ii] spectra that they saw in their [C ii] map as emission coming from the front and the backside of the PDR shell, but the optically thin [O i] 145 spectra show that the dip in the [C ii] spectra is due to self-absorption, not two separate velocity components. Neither can we distinguish more than one velocity component in the CO spectra, although there are areas, see Fig. 6 and Table 1, where the lines appear broadened, most likely due to contribution both from the front and the back wall of the PDR.

In most of the positions the [O i] 63 spectra show wings and line widths which are comparable to the [C ii] spectra, although they are often completely self-absorbed at the centre of the line. The [O i] 145 spectra are much narrower than the [O i] 63 |$\mu$|m and show profiles which are almost identical to the CO(7–6) emission, except at position −2″, −78″ where the CO(7–6) shows strong blue-shifted from the outflow powered by Sellgren’s star C (Sandell, private communication).

3.3 Column density of the atomic oxygen

The [O i] 63 |$\mu$|m emission is strongly self-absorbed over most of the reflection nebula, whereas the [O i] 145 is always observed in emission with no evidence for self-absorption. Even though [C ii] and [O i] 63 |$\mu$|m shows bright wings due to photoevaporation flows, the photoevaporating gas is not dense enough to be seen in [O i] 145 except in a few areas, see e.g. spectra at positions 40″, −40″, and 23″, −60″, where [O i] 145 shows a clear blue-shifted line wing. In order to derive the distribution of atomic oxygen based on the 63 and 145 |$\mu$|m spectra, we fit the [O i] 63 spectra using a two-slab toy model: the hot PDR shell (lying on either side of HD 37903 along the line of sight) is assumed to be the background layer, while the colder oxygen on the front side is taken to be the foreground layer. While the gas in this layer, where the PDR merges into the cold molecular cloud is largely molecular, most of the oxygen is still atomic. We further assume that the emission from the background layer is captured by the [O i] 145 emission. In this case the colder foreground oxygen primarily attenuates the background spectral emission at 63 |$\mu$|m by the factor exp [−τ₀exp [−4ln 2(v − v₀)²/Δv²]]. Here τ₀, v₀, and Δv denote the peak optical depth, the velocity at which τ₀ of the foreground cloud occurs and the full width at half maximum (FWHM) of the absorption profile due to the foreground material. As pointed out earlier, the [O i] 63 emission arising from a transition involving lower energy levels also trace the blue- and red-shifted gas in photoevaporation flows. The [O i] 145 emission, which is only excited in the warm and dense gas, does not show these broad-line wings. The background [O i] 63 emission is modelled by Gaussian profiles generated from fits to the observed [O i] 145 spectrum scaled by a factor of 2, which corresponds to the typical intensity ratio between the two [O i] lines in temperature units and is equal to I_63μm/I_145μm = 25 in energy units. Calculations for warm, medium-density PDRs with T_kin= 100 K and n=10⁵ cm⁻³ suggest a I_63μm/I_145μm ratio between 20–33 (fig. 10 Goldsmith 2019). Assumption of a higher ratio of 33 leads to an increase in the derived oxygen column density by 10–20 per cent of the value obtained in this work. Fig. 7 shows an example of the results of the fits to the [O i] 63 using the two-slab toy model at the selected positions and Table 2 presents the corresponding fitted parameters. At the selected positions, the central velocity of the absorbing foreground cloud lies between 10.1–10.7 km s⁻¹ which matches well with the velocity of the [O i] 145 emission spectra and have FWHM between 1.8–2.4 km s⁻¹. This suggests that the [O i] 63 absorption occurs in the colder foreground layers of the same PDR gas that accounts for the emission of [O i] 145. We have performed the fitting of the [O i] 63 profile over a larger number of positions in the map and obtained the distribution of the foreground cold absorbing component (Fig. 8).

The fits using the [O i] 145 spectra as the model for the background component do not reproduce the red- and blue-shifted wings observed in the [O i] 63 spectra, since the gas in the photoevaporation flows causing these wings is not dense enough to excite [O i] 145. We estimate the column density of the absorbing foreground layer by assuming that all the oxygen atoms reside in the lowest level of their ground state. The column density of all O atoms in the foreground layer can be then estimated from the centre opacity τ₀ and width (in km s⁻¹) of the [O i] 63 line by

The estimated column densities for the foreground gas ranges from 8 10¹⁷ to 2 × 10¹⁸ cm⁻². We can estimate a lower limit to the extinction if we assume that all the oxygen is atomic. With a typical gas abundance of 3.2 × 10⁻⁴ for oxygen and assuming most of the hydrogen to be molecular we estimate N(H₂) of the foreground gas to range between 1.2–3.2 10²¹ cm⁻², which translates to a minimum A_V of 1.3. However, since one third of the oxygen is tied up in CO in molecular gas a more realistic estimate for the A_V is 2 mag.

For the estimated column densities, based on non-LTE calculations using radex (Van der Tak et al. 2007), the [O i] 145 line is optically thin over a large range of temperatures (20–300 K) and volume densities (10⁴–10⁷ cm⁻³). This is also consistent with our observations. Since the T_mb ratio for the two [O i] lines depends on the physical conditions, the derived values of the peak optical depth are subject to the ratio assumed in this work. However, the fitted central velocity and line widths for the absorption profile are fairly robust against the assumed scale factor. Furthermore, in the estimate of N(O) we have assumed that all O atoms in the absorbing layer are in the ground ³P₂ level, which is reasonable because there is no hint of absorption in [O i] 145.

Since the [O i] 63 |$\mu$|m is optically thick, it cannot be used to estimate the column density of oxygen atoms in the background layer. We thus use the integrated intensities of the [O i] 145 |$\mu$|m spectral lines to estimate N(O). For this, we need an estimate of the temperature and density of the emitting gas. Since the [O i] 145 emission seems to follow the CO(11–10) emission quite closely, we can use the temperature and density estimated for CO(11–10), for which Sandell et al. (2015) estimated a kinetic temperature of 120 K and densities of 2 10⁵–1 10⁶ cm⁻³. Table 2 presents the column density of atomic oxygen at the seven selected positions that is required to produce the observed [O i] 145 emission (Table 1), estimated based on non-LTE calculations using radex (Van der Tak et al. 2007) for these kinetic temperatures and densities (Table 2). Fig. A1 shows a plot of the intensity of [O i] 145 line as a function of N(O) for T_k = 120 K and n(H₂) = (0.2–1) 10⁶ cm⁻³.

3.4 The morphology of the reflection nebula and the PDR illuminated by HD 37903

Sandell et al. (2015) proposed that the PDR illuminated by HD 37903 formed an expanding thin, hot shell which is roughly spheroidal with the south-eastern part being almost spherical because the expansion is slowed down by the dense surrounding molecular cloud and the shell being much more extended to the north-west. They interpreted the double peaked [C ii] spectra seen on the north-western side of the star as originating from the front and backside of the PDR shell. They argued that the reason for not seeing clear double split spectra on the south-eastern side, was because the expansion was slowed down by the dense molecular cloud. While their model appears reasonable, we now know that it has some serious flaws. Even though we do see double peaked [C ii] spectra north-west of HD 37903, see e.g. Fig. 5, [O i] 145 is single peaked and centred between the two [C ii] peaks. This confirms that what we see is self-absorption in [C ii], not two velocity components. The [C ii] lines are quite broad, so it is probably true that they are broadened by emission from both the front and the backside of the PDR, but they cannot be separated into two distinct velocity components. The foreground extinction is definitely lower in the north and north-west. We do see some blue-shifted [C ii] wings due to photoevaporation flows. This is consistent with the cloud emission being more red-shifted in both ¹³CO(3–2) and [C ii] indicating that it is behind the PDR. It appears more likely that the PDR is fan-shaped rather than ‘egg-shaped’ and that it is most likely open in the north and north-west allowing FUV photons to escape, hence creating a large reflection nebula dominating the emission north and north-west of HD 37903.

4 DISCUSSION

Analysis of the observed [O i] 63 and [O i] 145 spectra suggests that the strong absorption features seen in the [O i] 63 profiles are caused by colder atomic oxygen lying between the observer and the warm PDR emission excited by HD 37903. The [C ii] spectra are also affected by self-absorption, but it is far less severe and barely noticeable over most of the reflection nebula. At the peak (27″, 18″) of the [C ii] emission (Fig. 6), the self-absorption results in the shifting of the [C ii] spectra toward blue-shifted velocities compared to [O i] 145 and the high-J CO lines, which are unaffected by self-absorption. Even the [C ii] spectrum on HD 37903 appears to be somewhat affected by self-absorption.

The two-slab toy model, which we used in Section 3.3 to model the observed [O i] 63 absorption features, shows that the column density of atomic oxygen in the foreground absorbing gas is quite substantial, ∼10¹⁸ cm⁻². This is because most of the oxygen remains atomic way into the cloud, i.e. from the hot PDR layer up to about an A_V of 10, as already shown by Tielens & Hollenbach (1985), while most of the carbon becomes molecular, i.e. gets tied up in CO. The contribution to the oxygen column density from the surface of the molecular cloud, which only gets ionized by the general interstellar radiation that is a factor of 10⁴ or 10⁵ times weaker than the PDR illuminated by HD 37903, is negligible. The estimated [O i] column densities depend to some extent on the ‘model’ for the background emission, which have been constructed using the [O i] 145 spectra, but as we can see in Table 2 the computed column densities only depend weakly on the assumed gas densities.

The [O i] 63 spectral profiles seen in NGC 2023 are similar to those observed in other strongly self-absorbed regions like S1 in RhoOph (Mookerjea et al. 2021) and toward W3 (Goldsmith et al. 2021), where the absorption was estimated to be due to colder foreground gas with N(O) of 3 10¹⁸ and 2–7 10¹⁸ cm⁻², respectively. The method for estimating the background for the S1 PDR in Rho Oph was similar to what we used in this paper, while for W3, which was not observed in [O i] 145 it was instead estimated by fitting the [O i] 63 spectra for positions which showed almost no evidence for absorption. In the [O i] 63 map in NGC 2023 there are only a few positions which show a single emission peak without any trace of absorption, thus the approach to use the [O i] 145 emission which has been measured is more robust.

Since we detect the [¹³C ii] F = 2–1 line in essentially the entire area mapped in [C ii], it is possible to check whether the [C ii] emission is optically thick although the [C ii] map does not go deep enough to measure [C ii] to [¹³C ii] ratio in individual positions. We have therefore divided up the map into smaller regions of about 1 sq. arcmin, each containing about 50–70 spectra. When we did this we found that the ratio of intensities between [C ii] to [¹³C ii] within errors is ∼30–40 after correcting for the relative intensity of the F = 2–1 line, 0.625. Considering the isotope ratio ¹²C/¹³C of 70 for the Orion region, this suggests the [C ii] optical depth to be ∼1.5–2, consistent with the value of 2 derived by Sandell et al. (2015). We obtain an estimate of the minimum kinetic temperature T_kin(C⁺) from the Planck-corrected peak of the [C ii] spectrum assuming a beam filling factor of unity. The lower limit of T_kin(C⁺) thus estimated for the entire map ranges between 60–92 K. Assuming an optical depth of 2 and T_kin(C⁺), and the integrated [C ii] intensity to estimate N(C⁺) using the equation (26) from Goldsmith et al. (2012) modified for optically thick emission as follows:

where A_ul = 2.3 × 10⁻⁶ s⁻¹, T_kin is the gas kinetic temperature, the collision rate is C_ul = R_uln with R_ul being the collision rate coefficient with H₂ or H⁰, which depends on T_kin, and n is the volume density of H. Since the critical density of the [C ii] transition is n_cr=3000 cm⁻³, and since it is likely that most of the [C ii] detected could be at such densities along with some emission arising from clumps with densities exceeding 10⁵ cm⁻³, we assume a density of 10⁴ cm⁻³ to estimate the N(C⁺) density distribution of the region. For this calculation we use the kinetic temperatures of C⁺ estimated as above, a density of n=10⁴ cm⁻³, and excitation from C⁺–H₂ collisions, with R_ul = 3.8 × 10⁻¹⁰ cm³ s⁻¹. The value of N(C⁺) estimated in the region show a small range of values between (1.8–4.4)× 10¹⁸ cm⁻². We note that barring a few exceptions (Okada et al. 2015; Mookerjea et al. 2021) for the estimate of N(C⁺) often a kinetic temperature of 100–120 K is assumed in contrast to using an estimate from the observed [C ii] peak temperatures as we have done here. From equation (2), we find that for a temperature of 60 K the N(C⁺) estimated will be approximately twice the value that would be estimated for T_kin=120 K. As seen at the selected positions, throughout the entire region the ratio of column densities of C⁺ and O atoms (seen in emission) is consistent with the comparable to the solar [O]/[C] abundance ratio of 3.5 within the uncertainties (Table 2). This is also consistent with the structure of the PDR expected where the [O i] 145 emission arises from a dense hot region which also emits in [C ii], although [C ii] emission additionally arises from the more diffuse and colder PDR gas.

Based on the FIR dust continuum maps at 160 to 500 |$\mu$|m, we estimate the molecular hydrogen column density to be around 7 × 10²¹ cm⁻², while to the south, beyond the ridge and into the molecular cloud the values are typically ∼2 × 10²² cm⁻². The estimated H₂ and C⁺ column densities at the selected positions are consistent with the solar C/H abundance ratio of 3 × 10⁻⁴ (Table 2).

5 SUMMARY AND CONCLUSIONS

We have studied the geometry of the NGC 2023 reflection nebula using spectrally resolved observations of the 63 and 145 |$\mu$|m transitions of [O i] along with the [C ii] 158 |$\mu$|m transition which enable us to map the three-dimensional distribution of the photodissociated gas vis-a-vis the location of the illuminating source HD 37903. Combination of the velocity-resolved [O i] 145 spectra with the [C ii] spectra have allowed us to improve our understanding of the morphology of the region. We conclude that the PDR is fan-shaped rather than ‘egg-shaped’ and that it is most likely open in the north and north-west allowing FUV photons to escape, that creates a large reflection nebula dominating the emission north and north-west of HD 37903. The [C ii] emission detects PDR gas uniformly distributed with N(C⁺)∼0.9–|$2.2\, 10^{18}$| cm⁻² suggest the contribution of both diffuse and dense components, while the [O i] 145 emission detects dense PDR gas with N(O) between (2–10) 10¹⁸ cm⁻². Additionally, the self-absorbed profile of the [O i] 63 indicate the presence of cold low-excitation atomic oxygen with N(O)∼0.8–1.5 10¹⁸ cm⁻² associated with the nebula and lying in the foreground between the illuminating source and the observer that is most pronounced towards the southern and south-western part of the nebula. The self-absorbed [O i] 63 profiles observed in NGC 2023 together with similar results in other Galactic PDRs (Goldsmith et al. 2021; Mookerjea et al. 2021) suggest the need to exercise caution when modelling galactic and extragalactic [O i] 63 intensities, in particular if spectrally unresolved, without accompanying observations of the optically thin [O i] 145 line.

ACKNOWLEDGEMENTS

BM acknowledges the support of the Department of Atomic Energy, Government of India, under Project Identification no. RTI 4002. This research has made use of the VizieR catalogue access tool, CDS, Strasbourg, France (DOI : 10.26093/cds/vizier). The original description of the VizieR service was published in 2000, A&AS 143, 23″.

DATA AVAILABILITY

The new data observed with SOFIA that is presented here will be available from the SOFIA data archive (https://irsa.ipac.caltech.edu/applications/sofia/).

Footnotes

References

APPENDIX A: RESULT OF NON-LTE CALCULATIONS USING radex

Fig. A1 shows the result of the non-LTE radiative transfer model calculations to estimate N(O) assuming T_kin= 120 K and n(H₂) = 2 10⁵ and 10⁶ cm⁻³ performed using radex (Van der Tak et al. 2007).