C60 IN PHOTODISSOCIATION REGIONS

PABLO CASTELLANOS<fnref n="1" />, OLIVIER BERNÉ<fnref n="2" />,<fnref n="3" />, YARON SHEFFER<fnref n="4" />, MARK G. WOLFIRE<fnref n="4" />, AND ALEXANDER G.G.M. TIELENS<fnref n="1" />

Received 2013 August 26; accepted 2014 August 1; published 2014 September 25

ABSTRACT

Recent studies have confirmed the presence of buckminsterfullerene (C60) in different interstellar and circumstellar environments. However, several aspects regarding C60 in space are not yet well understood, such as the formation and excitation processes, and the connection between C60 and other carbonaceous compounds in the interstellar medium, in particular polycyclic aromatic hydrocarbons (PAHs). In this paper, we study several photodissociation regions (PDRs) where C60 and PAHs are detected and the local physical conditions are reasonably well constrained to provide observational insights into these questions. C60 is found to emit in PDRs where the dust is cool (Td = 20–40 K) and even in PDRs with cool stars. These results exclude the possibility for C60 to be locked in grains at thermal equilibrium in these environments. We observe that PAH and C60 emission are spatially uncorrelated and that C60 is present in PDRs where the physical conditions (in terms of radiation field and hydrogen density) allow for full dehydrogenation of PAHs, with the exception of Ced 201. We also find trends indicative of an increase in C60 abundance within individual PDRs, but these trends are not universal. These results support models where the dehydrogenation of carbonaceous species is the first step toward C60 formation. However, this is not the only parameter involved and C60 formation is likely affected by shocks and PDR age.

**Key words:** astrochemistry – infrared: ISM – ISM: lines and bands – ISM: molecules – photon-dominated region (PDR)

*Online-only material:* color figures

1. INTRODUCTION

Buckminsterfullerene, or C₆₀, was first discovered in laboratory experiments aiming to understand the spectroscopy of carbon chain molecules in the interstellar medium (ISM) and circumstellar envelopes (Kroto et al. 1985). Based on their study, the discoverers of C₆₀ concluded that the molecule corresponded to a new form of carbon organized as a truncated icosahedron, usually compared with the old style black and white “soccer” ball.

C₆₀ is a super-stable molecule and is considered the prototypical fullerene, cage-like molecules of pure carbon. Krätschmer et al. (1990) developed a method to produce bulk quantities of C₆₀ by evaporating graphite electrodes in a helium atmosphere. They also confirmed the structure of the molecule through X-ray diffraction and its infrared (IR) spectrum. The discovery of fullerenes has also opened the research of carbon nanotubes, another different form of carbon which combines the properties of graphite and fullerenes (Ijima 1991).

The values for the IR active modes are 526, 575, 1182, and 1429 cm⁻<fnref n="1" />, which correspond to 18.9, 17.4, 8.5, and 7.0 µm (Menéndez & Page 2000). While the values for the frequencies show a good agreement, reported intrinsic band strengths vary widely. Several theoretical and experimental works have been carried out (Chase et al. 1992; Fabian 1996; Choi et al. 2000; Iglesias-Groth et al. 2011) but the results have large differences.

Since its very discovery, C60 was considered a potential component for interstellar dust and speculated to be the carrier of some of the diffuse interstellar bands (DIBs). Its survival in ISM conditions is supported by its high stability. However, until recently, unequivocal detection has not been possible. A first tentative detection was the association of two, weak, far-red absorption bands with C<sub>60</sub>⁺ (Foing & Ehrenfreund 1994). This identification has been contested and the issues were never fully

resolved (Maier 1994). The main obstacle for the detection of fullerenes in emission stems from the fact that the mid-IR spectra of almost all interstellar sources are dominated by the vibrational spectrum of polycyclic aromatic hydrocarbons (PAHs). Because of this, in sources with strong PAH emission, small amounts of fullerenes are difficult to detect, with their emission hidden by the PAH bands and the continuum.

Cami et al. (2010) recognized C₆₀ and C₇₀ bands in the spectrum of the planetary nebula (PN) Tc 1, which has no strong PAH bands. After the detection of these transitions and their association with C₆₀ several more objects have been investigated, showing that C₆₀ exists in a variety of sources with different evolutionary stages and physical conditions. Most of these works have dealt with PNe, including galactic and extragalactic sources in the Magellanic Clouds (García-Hernández et al. 2011; Otsuka et al. 2012), proto-PNe (Zhang & Kwok 2011), and in circumstellar envelopes of binary post-asymptotic-giant-branch stars (Gielen et al. 2011). There has also been a detection toward a young stellar object (YSO) and a Herbig Ae/Be (Roberts et al. 2012), which represent isolated pre-stellar objects. C₆₀ was also detected in photodissociation regions (PDRs) associated with both reflection nebula (RN; Sellgren et al. 2010) and H II regions (Rubin et al. 2011).

The excitation mechanism of C₆₀ is not yet clear. This is an important question to address since it determines which physical conditions of its environment are traced by the bands. Two different mechanisms have been suggested: Cami et al. (2010) consider the band ratios from Tc 1 to be consistent with a thermal distribution of the excited states, deriving an excitation temperature of ~330 K. This mechanism requires that C₆₀ is not in the gas phase, but in the solid state or deposited on dust grains.

Another excitation mechanism, proposed by Sellgren et al. (2010), assumes that C₆₀ remains in the gas phase and the bands originate from IR fluorescence. This is widely accepted

<footnote n="1">¹ Leiden Observatory, P.O. Box 9513, 2300 RA Leiden, The Netherlands; pablo@strw.leidenuniv.nl</footnote>

<footnote n="2">² Université de Toulouse, UPS-OMP, IRAP, F-31400 Toulouse, France</footnote>

<footnote n="3">³ CNRS, IRAP, 9 Av. Colonel Roche, BP 44346, F-31028 Toulouse Cedex 4, France</footnote>

<footnote n="4">⁴ Department of Astronomy, University of Maryland, College Park, MD 20742, USA</footnote>
as the excitation mechanism for PAHs and consists in the absorption of a single UV photon which leaves the molecule highly excited and leads to a redistribution of the energy between the vibrational modes. While for NGC 7023 the reported 7.0/18.9 µm ratio is in good agreement with fluorescence models, a much lower value was reported for this ratio in the reflection nebula NGC 2023 (Sellgren et al. 2010). We note though that this ratio is very difficult to determine in spectra dominated by the PAH features as the underlying plateau emission is very strong and broad. This holds for both NGC 7023 and NGC 2023 and is compounded when the program PAHFIT is used to extract the intensity for highly blended, weak features as the adopted intrinsic Lorentzian profile (inappropriate for fluorescence from highly excited species; see Tielens 2008) locates much of the emission in ill-fitted wings. Moreover, thermal as well as fluorescence analyses of observed C₆₀ band ratios in space are hampered by the unknown intrinsic strength of these modes.

Recently, Berné et al. (2013) reported a detection of C<sup>+</sup><sub>60</sub> in NGC 7023 through four bands at 6.4, 7.1, 8.2, and 10.5 µm. This would support the idea that C<sub>60</sub> is in the gas-phase, at least in environmental conditions similar to those found in NGC 7023.

The formation of C<sub>60</sub> is also a subject of debate. It seems self-evident to consider that C<sub>60</sub> could be built-up from small (hydro)carbon chains in the C-rich ejecta of Tc 1 and other PNe whose spectra are dominated by C<sub>60</sub>. In H-poor environments (not the case of these PNe), simple formation paths analogous to the chemical routes described by Kroto et al. (1985) might lead to C<sub>60</sub> formation (Cherchneff et al. 2000). However, neither mechanism is efficient enough to account for the widespread detection of C<sub>60</sub> in the ISM.

Another hypothetical formation route in PNe is through photo-chemical processing of hydrogenated amorphous carbon (HAC; Bernard-Salas et al. 2012; Micelotta et al. 2012). The details of this mechanism require the presence of large HAC clusters (N<sub>C</sub> > 100) with a high H atomic fraction (i.e., not exposed to strong or continued UV radiation). The sudden exposure of these small compounds to UV photons is speculated to lead to dehydrogenation and aromatization. This would result in giant fullerene cages which will shrink to C<sub>60</sub> and C<sub>70</sub> through C<sub>2</sub> losses.

Finally, Berné & Tielens (2012) recognized that the abundance of C₆₀ increases rapidly near the illuminating star of the reflection nebula NGC 7023, while the abundance of PAHs decreases. They propose that UV processing of PAHs leads to fullerene formation. Under the influence of radiation, PAHs will first become dehydrogenated, leading to the formation of graphene sheets. A continued exposure to high energy photons is expected to remove carbon from within the structure, which will force the formation of pentagonal rings and the bending of the sheet. This process of destruction might lead to several different intermediate stages. Because of its high stability, C₆₀ is expected to be the photo-product with the longest lifetime, and therefore it can be present at ~0.01% of the elemental C near the star.

In this work, our goal is to contribute to the understanding of the origin and evolution of C₆₀ in environments where the presence of PAHs is also strong. We test the formation mechanism proposed by Berné & Tielens (2012) by comparing the variations of the abundance of C₆₀ and PAHs as a function of physical conditions (G₀, n<sub>H</sub>). We also consider the spatial variations of these abundances within the sources, in particular with respect to the position of strong sources of radiation.

<figure><img src="image_1.png" /><figcaption>**Figure 1.** Composite image of NGC 2023. In blue is the IRAC 8 µm image, showing the distribution of PAHs. In green is the PACS 70 µm map, which traces warm dust grains. In red is PACS 160 µm, tracing colder dust. The cross indicates the position of HD 37903 and the black box shows the IRS field that we will use for this work.</figcaption></figure>

(A color version of this figure is available in the online journal.)

We used data from IRS (Houck et al. 2004) on board the *Spitzer* Space Telescope (Werner et al. 2004) and the PACS instrument (Poglitsch et al. 2010) on the Herschel Space Observatory (Pilbratt et al. 2010), for three PDRs: NGC 2023 North, Ced 201, and RCW 49. We complement this analysis by including observations of the Horsehead nebula, IC 63, and Orion's Veil. For these latter sources we do not give a detailed analysis since there is no detection of C₆₀, with the exception of Orion's Veil which was studied by Boersma et al. (2012).

This paper is organized as follows. In Section 2, we give a description of the physical and chemical conditions of the sources, taken from the literature. In Section 3, we describe the data from IRS and PACS that we use and explain the data reduction process. In Section 4, we present the reduced spectra and feature maps obtained from IRS, as well as the derivation of the UV field intensity from the PACS for each source. We also describe the observed variations of C₆₀ and PAH abundances within each source. In Section 5, we discuss our results and compare them with the C₆₀ formation model of (Berné & Tielens 2012). Finally, in Section 6, we summarize the main results and conclusions of this work.

## 2. SOURCES

### 2.1. NGC 2023

NGC 2023 is at a distance of 330–385 pc and is illuminated by HD 37903. This star is usually classified as a B1.5 V star (Brown et al. 1994), although more recent works have concluded that it is a Herbig Ae/Be star of spectral class B2 Ve (Mookerjea et al. 2009) based on UV spectra and IR excess. This RN is part of the dark cloud L 1630 in Orion, which also hosts the famous Horsehead nebula (B33). Surrounding HD 37903 there is a small H II region with a size of 0.015 pc (Knapp et al. 1975; Pankonin & Walmsley 1976). A far-IR study by Harvey et al. (1980) showed that the molecular cloud L 1630 lies behind HD 37903. The age of the stars in NGC 2023 lies in the range 1–7 Myr with the most massive stars falling toward the lower range (López-García et al. 2013).

Figure 1 gives a general overview of NGC 2023. We can see that the shape of the reflection nebula is roughly circular,
but HD 37903 is not completely centered in it but displaced to the southeast. The radiation from HD 37903 is the main source of UV photons and creates a bubble-shaped cavity, particularly clear in the IRAC 8 µm image (Figure 1), which traces mostly PAH emission. Outside the cavity most of the PAHs become faint as expected because of the extinction of the UV field. Surrounding the cavity, several filaments or ridges can be observed, most clearly in PAH and warm dust emission. In particular, close to HD 37903, the southern ridge is the most luminous part of the nebula, hosting a high concentration of YSOs (Lada et al. 1991; Mookerjea et al. 2009). Other less luminous ridges can be seen to the north and the west of HD 37903. The black box in Figure 1 shows the IRS field that will be the focus of our analysis, which lies on top of the northern ridge of NGC 2023.

From observations of far-IR fine structure lines of [C II], [O I], and [Si II], Steiman-Cameron et al. (1997) have established that the environment in NGC 2023 is clumpy. Their observations of the northern PDR are not spatially coincident with the position of the IRS field we analyze, corresponding to a position closer to the HD 37903. For this position they fit a model with $G_0 = 2 \times 10^4$ and $n_0 = 3 \times 10^3$ cm$^{-3}$. Pilleri et al. (2012) give a model for the same IRS field we observe. They derive $n_0$ as a function of the distance to HD 37903 fitting a power law for the density profile, from the density in the diffuse region (a free parameter in their fit) up to a maximum density, which they take from Fuente et al. (1995), based on CN observations. The density ranges from $2.4 \times 10^3$ cm$^{-3}$ to $2 \times 10^4$ cm$^{-3}$ throughout the field.

Sellgren et al. (2010) reported the detection of three features of $C_{60}$ at 7.0, 17.4, and $18.9$ µm in NGC 2023. A recent work by Peeters et al. (2012) also detected extended $C_{60}$ emission through the $18.9$ µm feature in two locations inside the nebula, using the IRS instrument on board *Spitzer*. They found in their southern position that the observed spatial variations of $C_{60}$ and PAHs are consistent with Berné & Tielens (2012), with $C_{60}$ appearing closer to the position of HD 37903 than the PAH features. However, in their northern position, they find the opposite with $C_{60}$ peaking further away from the illuminating star than the PAHs, which they suggest may be due to geometrical effects.

## 2.2. Ced 201

Ced 201 is an unusual RN, since it is the result of a chance encounter of a molecular cloud with a runaway B9.5 V star (BD +69 1231). This was estimated considering a radial velocity difference between the cloud and the molecular cloud of ~12 km s$^{-1}$ (Witt et al. 1987). In most RNe, the illuminating star is formed within the cloud. Due to its velocity, the star in Ced 201 not only affects the molecular cloud through the radiation field, but also induces shocks and turbulence. An arc-like structure can be observed to the east of the star, which is particularly clear in white in Figure 2. This structure has been interpreted as a shock due to the velocity difference of the star and the molecular cloud. However, no signature of a shock has been observed in the line profiles of CO (2–1) (Cesarsky et al. 2000). The molecular cloud is part of the Bok Globule B175 (L 1219) at a distance of ~400 pc (Casey 1991).

Figure 2 shows the region of Ced 201. The star is entering the molecular cloud from the west, generating the aforementioned arc-like structure. We can see that the large dust grains (in red) and with them the molecular cloud itself, are located to the top of the image. The hot dust traced by MIPS 24 µm is located

<figure><img src="image_2.png" /><figcaption>**Figure 2.** Composite image of Ced 201. Blue corresponds to the IRAC 8 µm image, green to MIPS 24 µm, and red to PACS 70 µm. The box is the IRS field that will be used in this work. The cross marks the position of the illuminating star BD +69 1231.</figcaption></figure>

(A color version of this figure is available in the online journal.)

closer to BD +69 1231. PAH emission, traced by IRAC 8 µm in Figure 2 (blue), is seen mostly located near the star and to the south. PAHs and large dust grains coexist in the proximity of BD +69 1231.

Given that the PDR is the product of a chance encounter, we can estimate its age by the measured velocity of the star. In order to do this, we use the measured proper motion of the star (Høg et al. 2000) and we calculate the angular distance from the star to the edge of the cloud in the direction of the proper motion. This gives us a rough estimate of ~1500 yr.

Young Owl et al. (2002) used far-IR fine structure line intensities of [C II] and [O I] (for the latter only upper limits were derived) to determine that the gas physical conditions are $n = 4 \times 10^2$ cm$^{-3}$ and $T = 200$ K, while for the radiation field they find $G_0 \sim 300$. On the other hand, a previous study by Kemper et al. (1999) found higher values of $n = 1.2 \times 10^4$ cm$^{-3}$ and $T \sim 330$ K. They also find $n(H_2) = (5 \pm 1) \times 10^3$ cm$^{-3}$. These last values are derived from a PDR model that fits the observed fine structure atomic lines (again [C II] and upper limits of [O I]) and molecular sub-millimeter lines. They detected CO, $^{13}$CO, and HCO$^+$ and in addition put upper limits for CS and $C^{18}$O. The discrepancy in the values for the density stems from the fact that that Kemper et al. (1999) based their model not only on the fine-structure lines but also on molecular data. We expect that their values represent a better constraint on the density values. This is also supported by Casey (1991), who finds $n > 2000$ cm$^{-3}$ based on dust emission from far-IR.

## 2.3. RCW 49

RCW 49 is one of the most luminous and massive H II regions in the galaxy. The region is powered by the compact cluster Westerlund 2 (Westerlund 1960) where several OB stars have been detected, including at least twelve O-stars earlier than O7, the earliest of type O3 V(f) (Rauw et al. 2007). It also hosts two Wolf–Rayet stars, WR20a and WR20b. WR20a was found to be a binary system of two WN6ha stars with individual masses of ~70 $M_{\odot}$ (Rauw et al. 2004). WR20b appears to be a single star with spectral type WN6ha. The presence of Wolf–Rayet stars imply a cluster age of 2–3 Myr (Piatti et al. 1998). The molecular density and kinetic temperature have been determined through $^{12}$CO and $^{13}$CO observations by Ohama et al. (2010), with the
<figure><img src="image_3.png" /><figcaption>**Figure 3.** Composite image of RCW 49. Blue corresponds to the IRAC 8 µm image, while green and red correspond the PACS images at 70 µm and 160 µm, respectively. The boxes show the seven IRS fields which will be used in this work. The numbers on top of each box correspond to the nomenclature used from here on to refer to the individual regions. The crosses indicate the position of both Wolf-Rayet stars, with WR20a to the north and WR20b to the south. Finally, the dashed circle indicates the position of the cluster Westerlund 2, which hosts most of the OB stars in the region.</figcaption></figure>

(A color version of this figure is available in the online journal.)

temperature ranging from 30 to 150 K and a typical density of $n_{H_2} \sim 3000 \text{ cm}^{-3}$. Y. Sheffer et al. (in preparation) estimate the density in five of the IRS fields that we consider in this work. They fit models to the H$_2$ IR line intensities in order to get the density parameters. These values will be discussed in detail in Section 4.4.

Figure 3 shows the environment of RCW 49. PAHs and dust are seen all over the nebula, with the brightest ridges near the Wolf-Rayet stars. The overall shape is somewhat elongated in the north-south direction, resembling a kidney. Two prominent bubbles are observed in radio continuum images, one surrounding Westerlund 2 and the other around WR20b. The bubble around Westerlund 2 is open to the west (Whiteoak & Uchida 1997). In this region, dust and PAHs coexist with the ionized gas or are at least embedded in neutral gas mixed with the ionized region (Churchwell et al. 2004). PAH emission traced by IRAC 8 µm (blue in Figure 3) is more intense in the surroundings of Westerlund 2. The ridges, bubbles, and pillars that can be seen in Figure 3 are evidence of the strong interaction of the parental molecular cloud with the stellar radiation and winds from Westerlund 2 and the more recently formed stars. As we see, the seven IRS fields considered in this work cover regions with very different conditions, allowing us to probe a large range of UV field intensities. In the following, when we refer to individual positions we will use the numbers shown above the respective box in Figure 3.

The distance to RCW 49 has been a subject of debate, with different values ranging from 2.3 to 7.9 kpc (Brand & Blitz 1993; Whiteoak & Uchida 1997; Moffat et al. 1991), where the lower limit comes from radio continuum observations and kinematic studies and the upper limit is derived from the luminosity distance. A recent photometric analysis by Carraro et al. (2013) derives the reddening and extinction law toward several members of Westerlund 2. Using this, they correct the apparent distance modulus in order to obtain the distance to the cluster, $d = 2.85 \pm 0.43$ kpc.

## 2.4. Additional Sources

We consider additional observations of Orion's Veil (Boersma et al. 2012) and NGC 7023 (Berné & Tielens 2012), for which $C_{60}$ has been detected. In the case of Orion's Veil, there are 11 positions available, which are labeled with increasing distance to Orion's Bar as I4–1, M1–4, and V1–3, respectively. From these we do not use the two closest regions, I4 and I3, given that in the corresponding PACS images they suffer from strong contamination with the emission from the Bar. We also discard the farthest region, V3, considering that it looks edge on to the back PDR, and thus does not probe exactly the same environment as the remaining regions (Boersma et al. 2012).

We also considered for our results and discussion two PDRs for which the $C_{60}$ is not detected: IC 63 and the Horsehead nebula. In the case of the Horsehead, we used data corresponding to the “mane,” where there is no detection of any excess at about 19 µm (where the strongest $C_{60}$ feature falls) over the continuum, which is strong in this region. We consider as upper limits the highest value that could be fitted by PAHFIT without creating a notorious bump (always fitted as a detection below $3\sigma$). For IC 63, there are areas where emission over the continuum that cannot be fitted by the [SIII] line alone. These pixels are fitted by PAHFIT with confidence over $3\sigma$. Upon visual inspection, however, we find that the fit is of poor quality and we consider that $C_{60}$ is not detected there. Because of this, we consider the fit as an upper limit only. Physical conditions for these additional sources will be provided in Sections 4.3 and 4.4.

## 3. DATA

### 3.1. IRS Data

We used data from the Infrared Spectrograph (IRS; Houck et al. 2004) on board the Spitzer Space Telescope (Werner et al. 2004). The spectrograph consists of four modules according to the wavelengths and resolution: short wavelengths with low resolution (SL), long wavelengths with low resolution (LL), short wavelengths with high resolution (SH), and long wavelengths with high resolution (LH).

For the low-resolution modules the wavelength coverage goes from 5.2 µm to 38.0 µm for the combination of SL and LL and have a spectral resolution, $R = \lambda/\Delta\lambda$, between ~60 and ~130. In the case of the high-resolution modules, the combined coverage of SH and LH ranges from 9.9 µm up to 37.2 µm with a spectral resolution $R \sim 600$.

We gathered the spectral cubes from the Spitzer Heritage Archive.<fnref n="5" /> For NGC 2023 and Ced 201 we used observations from the program with PID 3512 (PI C. Joblin). In the case of NGC 2023, we considered observations performed with SL and LL modules with AOR keys 12014848 and 12011264, respectively. For Ced 201 we selected again SL and LL modules, in this case with AOR keys 11047936 and 11047680, respectively.

For RCW 49 we have seven fields, all of them in SH module, from the program with PID 20012 (PI M. Wolfire). The AOR keys of each of them are 13812992, 13813248, 13813504, 13813760, 13814016, 13814272, and 13814528. These regions are not contiguous and will be referred to according to the numbers shown in Figure 3.

The data reduction and cube generation were performed using the CUBISM software (Smith et al. 2007a). The recommended data reduction begins by taking a background when available.

<footnote n="5">⁵ http://sha.ipac.caltech.edu</footnote>
<figure><img src="image_4.png" /><figcaption>**Figure 4.** Examples of a fitted low-resolution spectrum of NGC 2023. The black line with error bars is the observed spectrum while the red line is the total fitted spectrum. The individual components are shown in blue (PAH and C₆₀ features), green (unresolved H₂ lines), and orange (modified black body continua).
(A color version of this figure is available in the online journal.)</figcaption></figure>

None of the observations has a dedicated background observation. The low-resolution modes have “outrigger” slits that fall outside the area of interest. In the case of Ced 201, this slit falls on an empty patch of sky so we used this to do the background subtraction. In the case of NGC 2023, the slit falls in regions with important emission from the same source which leads us not to consider a background for this source. In this case and in RCW 49, the lack of background should not be a very important issue since both sources are very bright.

## 3.2. PACS Data

In order to measure the FIR intensities, we used photometric data from the Photodetector Array Camera & Spectrometer (PACS; Poglitsch et al. 2010) on board the ESA Herschel Space Observatory (Pilbratt et al. 2010). In all three sources we retrieved the available data at 70 and 160 µm from the Herschel Science Archive.<fnref n="6" /> The observations of NGC 2023, Ced 201, and RCW 49 have respective OBSIDs of 1342227049/50, 1342196809/10, and 1342255009/10.

# 4. ANALYSIS

## 4.1. Spectral Fit

In order to fit the different features from the IRS spectra, we used the tool PAHFIT (Smith et al. 2007b). This tool fits a stellar continuum and several modified black bodies at different temperatures. These black bodies are used to fit the mid-IR continuum, but their temperatures do not represent the temperature of the dust, since grains emitting at these wavelengths do not reach an equilibrium with the radiation field, but are rather transiently heated to high temperatures (Draine 2003). It also considers several unresolved atomic lines as well as H₂ rotational lines. Finally it fits different PAH features and the 18.9 µm feature.

In Figure 4, we present an example of a low-resolution spectrum from NGC 2023 from the area where the 18.9 µm feature peaks. The fit matches reasonably well the observed spectrum with the exception of the PAH bands between 13 and 14 µm, for which the peak positions seem to be somewhat different from the values considered by PAHFIT. However, these bands are rather weak and when considering the integrated PAH band intensities their contribution will be small. These

<figure><img src="image_5.png" /><figcaption>**Figure 5.** Example of a fitted low-resolution spectrum from Ced 201. The black line with error bars is the observed spectrum while the red line is the total fitted spectrum. The individual components are shown in blue (PAH and C₆₀ features), green (unresolved H₂ lines), and orange (modified black body continua).
(A color version of this figure is available in the online journal.)</figcaption></figure>

characteristics are also observed in the spectrum of Ced 201 (Figure 5). In both sources, the 18.9 µm feature associated with C₆₀ is clearly present. The only unresolved lines fitted in these spectra are from H₂ pure rotational lines.

The main PAH features at 6.2, 7.7, 8.6, 11.2, and 12.7 and the plateau at 17.0 µm are also prominent. This allows us to consider in our analysis a pixel by pixel comparison of the variations of C₆₀ and PAHs. There are however some residual problems for the PAH bands: for instance, in NGC 2023 there is a dip around 10 µm which may due to saturation of the peak-up array as described in the IRS handbook<fnref n="7" /> (Section 7.3.5). In both sources, the issues are not relevant for the analysis since we will be focusing on the individual features as fitted by PAHFIT and no such feature falls in the wavelengths with problems.

In the case of the 18.9 µm feature, we found that in some pixels, the total intensity is underestimated. This defect is observed in NGC 2023, for pixels that lie to the south of the field. The origin of this problem comes from the continuum, which has a drop after 20 µm and thus is not properly fitted, affecting in turn the fit of the 18.9 µm feature. In order to circumvent this issue and check the accuracy of the fit in other pixels, we consider a local linear continuum around 18.9 µm and later fit a Drude profile for the feature. This procedure confirms that the only significant differences between the result of PAHFIT and our fit for the 18.9 µm feature appear in the aforementioned region. For these pixels, we use the result from our local fit instead of PAHFIT. Other features in the area are not significantly affected by the misplaced continuum.

One difference between the two sources is the feature-to-continuum ratio. The PAH features in Ced 201 have a significantly smaller intensity ratio with respect to the continuum than what we observe in NGC 2023. This is not the case for the 18.9 µm C₆₀ feature, which has a similar feature-to-continuum ratio in both sources. While hard to disentangle from the H₂ S(1) line and the 17.0 µm plateau, the 17.4 µm band also seems to have similar feature-to-continuum ratio in both sources. This could indicate a large contribution of C₆₀ to this interstellar band (Sellgren et al. 2010; Berné & Tielens 2012; Peeters et al. 2012).

The high-resolution spectra were also analyzed using PAHFIT. However, we have to take into account that PAHFIT was created and optimized for fitting low-resolution spectra, such as the ones from Ced 201 and NGC 2023. In the case of RCW 49, the spectra

<footnote n="6">6 http://herschel.esac.esa.int/Science_Archive.shtml</footnote>

<footnote n="7">7 http://irsa.ipac.caltech.edu/data/SPITZER/docs/irs</footnote>
<figure><img src="image_6.png" /><figcaption>**Figure 6.** Example of fitted a high-resolution spectrum from RCW 49, Region 2. The black line with error bars is the observed spectrum while the red line is the total fitted spectrum. The individual components are shown in blue (PAH and C₆₀ features), green (unresolved H₂ and atomic lines), and orange (modified black body continua).</figcaption></figure>

(A color version of this figure is available in the online journal.)

were taken in SH mode for all the regions considered here. This forces us to change the input parameters in PAHFIT. Moreover, for the unresolved atomic and molecular lines, we changed the width of the Gaussians to match the improved resolution of the SH mode. We also modified some of the PAH features to improve the fit.

The fitting procedure follows the same principle as in the low-resolution mode. Besides the previously mentioned modifications of the parameters, we were also forced to remove all features that fall outside the range covered by the SH mode since this causes problems with the fitting routine. If this last modification is not implemented, then the wings of the adopted Drude profiles of the out-of-range PAH features are sometimes fitted as part of the continuum, giving wrong fits.

We also found that this modified fit runs into problems when the feature to continuum ratio becomes small. This is the case for region 1 in RCW 49. In particular, at longer wavelengths the fit tends to either over- or underestimate the continuum, resulting in poor fits for the 18.9 µm feature. Furthermore, the 18.9 µm feature is not visually recognizable in the spectra. However, the continuum at this wavelength is very high for this region. Even the presence of an 18.9 µm feature with an intensity comparable with the highest values detected among our sources would not be recognizable. For these reasons we decided to exclude region 1 from our analysis.

We searched for potential errors in the fits for the continuum but in all cases the fit is consistent with the data. In Figure 6, we present an example of region 2 in RCW 49, which has the strongest 18.9 µm feature. The 18.9 µm feature varies in intensity between the different regions, but can be recognized in the spectra in most pixels. In some cases, the detection is only marginal, but in our final analysis we exclude all points with detection below 3σ with respect to the error bars given by PAHFIT.

## 4.2. Spectral Maps

Once the spectrum for each pixel have been fitted, we build maps for each of the PAH features, the spectral lines and the 18.9 µm feature. We also create maps for the errors, as derived by PAHFIT. These maps give the integrated intensity for each spectral feature and will be the basis for our analysis. In order

to restrict our analysis to those points with clear detection, we also flagged pixels for which the fitted intensity falls below a 3σ noise level. Finally, in order to get smoother images and correct for cases in which particular pixels have large differences in flux with respect to its neighbors, we performed a median filter with a 2 × 2 kernel.

In Figure 7, we present the maps of the 11.2 µm PAH feature, the 18.9 µm feature corresponding to C₆₀ and the H₂ 0–0 S(1) line at 17.0 µm. We show the results for NGC 2023, Ced 201, and Reg. 4 of RCW 49. We can see that in all the cases the PAH feature peaks closer to the illuminating star than the C₆₀ or H₂ lines. This contrasts with the observations on NGC 7023, where C₆₀ peaks closer to the illuminating star than the PAHs.

In NGC 2023, we identify a bar-like structure in all the maps, running from the southeast to the northwest. However, the exact position of the bar varies: on the 11.2 µm map it is displaced to the south when compared with the 18.9 µm or H₂ S(1) maps. Other variations observed in the bar concern the position of the peak intensity, with 18.9 µm displaying two peaks located toward the northwest and southeast of the bar, while H₂ S(1) shows one peak located roughly in the middle of the bar and 11.2 µm showing also one peak, but displaced to the northwest.

In the case of Ced 201, there is a clear asymmetry in the cases of the H₂ and C₆₀ maps. This asymmetry is also present but in a much more contained way in the case of the PAH maps. In the case of H₂ and C₆₀, we can see in both maps that the arc toward the northwest appears much more clearly than in the PAH maps. Also in both cases, the peak intensity appears to the northeast of the star, in the position of the main part of the molecular cloud as shown in Figure 2. In contrast, the PAH features are centered near the star, and, for the 8.6 and 7.4 µm features, the coincidence with IRAC 8 µm is very clear.

For RCW 49, we can see that the H₂ peaks normally avoid the peaks of the PAH features, and are typically found further away from the illuminating star. For the 18.9 µm feature, we can see that, in most regions we have little variations in intensity. As is the case for NGC 2023 and Ced 201, 18.9 µm has a clearly different spatial distribution when compared to the PAH features in all the regions. In RCW 49, all the PAH features follow a similar pattern, with minor variations in peak positions and relative strengths. However, when comparing to the 18.9 µm maps, we find much stronger variations in the peak position and relative intensities in most of the sources. Sometimes the PAH and C₆₀ morphology show a superficial resemblance, but the relative strength of the different components is very different.

In summary, for all the sources where we detect the 18.9 µm feature, its spatial distribution noticeably different with respect to the distribution of the PAH features. Although the different PAH bands show spatial variations with respect to each other, they resemble each other much more closely than the 18.9 µm feature. In most cases, the peak intensity of the C₆₀ band is seen farther from the central star than the PAH bands. This behavior is the opposite of what is observed in NGC 7023 (Berné & Tielens 2012) and, in principle, would not support the formation of C₆₀ by UV processed PAHs. It has been suggested that the formation of C₆₀ in the ISM starts with the dehydrogenation of PAHs which is a balance between UV photolysis and reactions with H, we will consider in the next two subsections the intensity of the local radiation field and the hydrogen density.

## 4.3. Calculation of G₀

Considering that PAH and C₆₀ are stochastically heated by UV-photons, we will need the intensity of the UV field in order
<figure><img src="image_7.png" /><figcaption>**Figure 7.** Intensity maps of the selected features: H₂ 0-0 S(1) line, 11.2 µm PAH feature, and 18.9 µm C₆₀ feature, all convolved to the PACS 70 µm beam. The top panels correspond to NGC 2023, the middle ones to Ced 201, and the bottom panels to Reg. 4 in RCW 49. The color scale is normalized with respect to the maximum intensity in the corresponding frame. The offsets are with respect to the position of the illuminating star, toward which the arrow or cross also points (for RCW 49 we chose the position of WR20a).</figcaption></figure>

(A color version of this figure is available in the online journal.)

to have a proper measure of the abundances. The intensity of the mid-IR bands will be proportional to the product of the column density of PAHs or C₆₀ with the UV field intensity.

We will determine the FUV field intensity from the FIR spectrum of dust in terms of G₀, which is a measure of the FUV field in terms of the Habing field (1.6 × 10⁻³ erg cm⁻² s⁻¹; Habing 1968). It considers the flux of photons with energies between 6 eV and 13.6 eV. One way of calculating G₀ is to consider that in a PDR nearly all the radiation is absorbed by dust in the UV and re-emitted in the far infrared (FIR). Assuming this and that grains act like a modified black body (Tielens 2005, Section 9.7), we can write

$$G_0 = 8.3 \times 10^3 \frac{V}{\ell S} \tau_{v_0} v_0^{-\beta} \int v^{\beta} B(v, T_d) dv, \quad (1)$$

where v₀ is a reference frequency where the optical depth is τv₀, β is the spectral index, and T_d is the dust temperature. The factor V/ℓS accounts for the geometry of the cloud, S is the surface area that is facing the illuminating star, V is the volume of the cloud, and ℓ corresponds to the length along the line of sight. The problem with this method is that the three-dimensional geometry is generally not well known, introducing uncertainties. For example, for an edge-on disk the geometric factor would be equal to the ratio of the disk thickness over ℓ while for a sphere this factor is equal to one. Since this method relies on the emission from dust, it will underestimate the value of G₀ in regions devoid of dust. In the three PDRs considered here, the environment is dusty and the IRS fields are away from hard sources of radiation, so we expect our G₀ estimate to provide a reasonable approximation.
<figure><img src="image_8.png" /><figcaption>**Figure 8.** Map of $G_0$ calculated from the FIR intensity for NGC 2023 with an assumed spherical geometry, thus setting $V/\ell S = 1$. The cross shows the position of the illuminating star and the square the position of the IRS field.
(A color version of this figure is available in the online journal.)</figcaption></figure>

To calculate the total FIR intensity we use the PACS images at 70 and 160 µm from each source. In evaluating Equation (1), we have corrected for the emission missing in the 70 and 160 µm filters by fitting modified black-bodies to those data. For this, we first matched the 70 µm images to the 160 µm resolution. With this we then fitted a modified black-body to these two points, on a pixel by pixel basis. The values to be fitted are the dust temperature $T_d$ and the optical depth $\tau_{v0}$. We fixed $v_0$ at 1000 GHz and took $\beta = 1.7$. The value of $\beta$ in the ISM is usually found to be somewhere between one and two, with some dependency on the dust temperature (Dupac et al. 2003). It has been observed that for $T_d \lesssim 20$ K then $\beta \sim 2$, while for $T_d \gtrsim 30$ K a value of $\beta \sim 1.5$ is found (Shetty et al. 2009). We chose an intermediate value since in all our sources we expect dust with temperatures in both regimes. We are not fitting $\beta$ directly since we only have two points for our fit which is not enough to realistically fit three parameters. To make the fit we used a Levenberg-Marquardt fitting routine. With this fit, we derive the values of $T_d$ and $\tau_{v0}$ which we can later use to determine $G_0$.

We assumed spherical geometry in the three cases. For NGC 2023 (Figure 8) this is likely a good assumption since the RN is almost circular, even though it is not exactly centered on the illuminating star. An additional factor must be added in Equation (1) to account for the fact that, particularly for B-type stars, there is a significant part of the radiation at wavelengths longer than 6 eV that will heat the dust but does not correspond to the FUV as defined by $G_0$. Based on the spectral type, we use a factor of 0.7 for the fraction of the stellar luminosity falling in the range between 6 and 13.6 eV (Steiman-Cameron et al. 1997). We check our calculation against the result of Sheffer et al. (2011). They found $G_0 = 1.7 \times 10^4$ in NGC 2023 S, while the region we are focusing on is about 2.3 times farther out. This is equivalent to a geometric dilution factor of 5.3, giving an expected value of $G_0 = 3200$, which is in agreement with the values seen in Figure 8.

In the case of Ced 201 (Figure 9), the spherical cloud assumption is not entirely realistic since we have a chance encounter with the molecular cloud and the star is likely situated at the edge of the cloud. We can see that in this case the peak of $G_0$ is displaced from the illuminating star and, furthermore,

<figure><img src="image_9.png" /><figcaption>**Figure 9.** Map of $G_0$ calculated from the FIR intensity for Ced 201 assuming spherical geometry ($V/\ell S = 1$). The cross shows the position of the illuminating star and the square corresponds to the IRS field.
(A color version of this figure is available in the online journal.)</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>**Figure 10.** Map of $G_0$ calculated from the FIR intensity for RCW 49 assuming spherical geometry ($V/\ell S = 1$). The crosses show the position of the illuminating stars and the squares the position of the five IRS fields for which we have density estimates.
(A color version of this figure is available in the online journal.)</figcaption></figure>

the region with highest $G_0$ corresponds to the arc structure seen northeast from BD +69 1231 in the molecular cloud, with much less intensity to the southwest. In this case, the fraction of the stellar luminosity emitted between 6 and 13.6 eV corresponds to 0.3 (Young Owl et al. 2002).

In the case of RCW 49 (Figure 10), the exact three-dimensional geometry is not well known, but it looks approximately circular and we expect that the spherical approximation will not introduce large errors. We will take into account in our discussion the possible effects that arise from changing the geometrical factor. Since the illuminating source of RCW 49 is dominated by the O-type stars we consider a factor of one for the fraction of the luminosity emitted in the FUV. In reality, a significant fraction will be emitted at wavelengths shorter than 13.6 eV in this case, but these photons will be absorbed and downgraded in the HII region.
<table><caption>**Table 1**
Average Physical Conditions and C<sub>60</sub> Fractions</caption><thead><tr><th>Region</th><th>n<sub>H</sub><br>(cm<sup>-3</sup>)</th><th>G<sub>0</sub></th><th>G<sub>0</sub>/n<sub>H</sub></th><th>f<sub>C,C<sub>60</sub></sub><br>(10<sup>-5</sup>)</th><th>Ref.</th></tr></thead><tbody><tr><td>NGC 2023</td><td>2.4 &times; 10<sup>3</sup></td><td>5200 &plusmn; 200</td><td>2.1</td><td>4.2</td><td>1</td></tr><tr><td>Ced 201</td><td>1.2 &times; 10<sup>4</sup></td><td>600 &plusmn; 60</td><td>0.050</td><td>44</td><td>2</td></tr><tr><td>RCW 49-2</td><td>4.0 &times; 10<sup>3</sup></td><td>11000 &plusmn; 1300</td><td>2.7</td><td>13</td><td>3</td></tr><tr><td>RCW 49-4</td><td>1.0 &times; 10<sup>4</sup></td><td>20000 &plusmn; 3000</td><td>2.0</td><td>7.2</td><td>3</td></tr><tr><td>RCW 49-5</td><td>2.5 &times; 10<sup>4</sup></td><td>3000 &plusmn; 400</td><td>0.12</td><td>11</td><td>3</td></tr><tr><td>RCW 49-6</td><td>2.5 &times; 10<sup>4</sup></td><td>5100 &plusmn; 300</td><td>0.21</td><td>12</td><td>3</td></tr><tr><td>RCW 49-7</td><td>2.5 &times; 10<sup>4</sup></td><td>34000 &plusmn; 2500</td><td>1.4</td><td>5.5</td><td>3</td></tr><tr><td>NGC 7023</td><td>1.5 &times; 10<sup>2</sup></td><td>10000 &plusmn; 7000</td><td>60</td><td>9.9</td><td>4</td></tr><tr><td>Veil-I2</td><td>8.7 &times; 10<sup>3</sup></td><td>31000 &plusmn; 1900</td><td>3.6</td><td>56</td><td>5.6</td></tr><tr><td>Veil-I1</td><td>8.1 &times; 10<sup>3</sup></td><td>22000 &plusmn; 500</td><td>2.7</td><td>41</td><td>5.6</td></tr><tr><td>Veil-M1</td><td>5.1 &times; 10<sup>3</sup></td><td>15000 &plusmn; 600</td><td>2.9</td><td>36</td><td>5.6</td></tr><tr><td>Veil-M2</td><td>5.1 &times; 10<sup>3</sup></td><td>13000 &plusmn; 600</td><td>2.5</td><td>32</td><td>5.6</td></tr><tr><td>Veil-M3</td><td>5.1 &times; 10<sup>3</sup></td><td>8700 &plusmn; 400</td><td>1.7</td><td>29</td><td>5.6</td></tr><tr><td>Veil-M4</td><td>3.4 &times; 10<sup>3</sup></td><td>5400 &plusmn; 200</td><td>1.6</td><td>30</td><td>5.6</td></tr><tr><td>Veil-V1</td><td>1.7 &times; 10<sup>3</sup></td><td>3500 &plusmn; 100</td><td>2.0</td><td>35</td><td>5.6</td></tr><tr><td>Veil-V2</td><td>2.6 &times; 10<sup>3</sup></td><td>5200 &plusmn; 100</td><td>2.0</td><td>31</td><td>5.6</td></tr><tr><td>Horsehead</td><td>8.0 &times; 10<sup>3</sup></td><td>~100</td><td>0.013</td><td>&lt;4.3</td><td>7</td></tr><tr><td>IC 63</td><td>3.0 &times; 10<sup>3</sup></td><td>~650</td><td>0.22</td><td>&lt;6.0</td><td>8.9</td></tr></tbody></table>

**Notes.** For NGC 2023, Ced 201, RCW 49, and the Veil, the values of G<sub>0</sub> were calculated following the procedure described in Section 4.3. Other values of G<sub>0</sub> and n<sub>H</sub> are taken from the given references. For NGC 7023 we use lower limits for G<sub>0</sub> and present the values at the peak of C<sub>60</sub> abundance.

**References.** (1) Pilleri et al. (2012); (2) Kemper et al. (1999); (3) Y. Sheffer et al. (in preparation); (4) Berné & Tielens (2012); (5) Boersma et al. (2012); (6) Rubin et al. (2011); (7) Abergel (2003); (8) Jansen et al. (1994); (9) Thi et al. (2009).

For the additional PDRs included in this work we only used the far-IR calculation of G<sub>0</sub> for Orion's Veil, resulting in G<sub>0</sub> ranging from ~30000 for the I2 field (the closest to the Bar) to G<sub>0</sub> ~ 3000-5000 for the fields that lie farther from the Bar. For NGC 7023 we used the values derived by Berné & Tielens (2012), which result in a lower estimate of G<sub>0</sub> = 10000 ± 7000 at the peak of the 18.9 µm feature. The Horsehead nebula was studied by Abergel (2003), who find G<sub>0</sub> ~ 100. For IC 63 the estimate by Jansen et al. (1994) is of G<sub>0</sub> ~ 650.

Finally, comparing Figures 1-3 against Figures 8-10 we see that the general structure seen in the composite images is clearly apparent in the G<sub>0</sub> distribution. The G<sub>0</sub> maps show the most characteristic features in the case of each PDR, such as the southern and northern ridges in NGC 2023, the bright bar observed in RCW 49 between the positions of WR20a and WR20b, and the arc-like structure seen in Ced 201 to the northeast of BD +69 1231. The derived dust temperatures inside the IRS fields of the three PDRs range between 20 and 40 K, which is much too low for emission from C<sub>60</sub> embedded in grains.

### 4.4. The Atomic H Density

In order to estimate the variations in abundance of PAHs and C<sub>60</sub> with respect to the physical conditions of each source, we need to derive the atomic hydrogen densities, n<sub>H</sub>, for each region. We estimate these based on previous studies and we derive a single value as representative for each of the regions we consider. A summary of the results for each field is presented in Table 1.

For NGC 2023 we base our estimate on the densities calculated by Pilleri et al. (2012). Their model uses a varying density inside the IRS field, which ranges from 2.4 × 10<sup>3</sup> cm<sup>-3</sup> to

2 × 10<sup>4</sup> cm<sup>-3</sup>. We will consider their value for the density in the inner part of the field, which will dominate the emission.

In the case of Ced 201, we base our estimate of n<sub>H</sub> on the values for total hydrogen density and molecular densities from Kemper et al. (1999) and Young Owl et al. (2002). The G<sub>0</sub> value derived from the IR measurement is in good agreement with the value derived from the PDR/molecular cloud analysis of Kemper et al. (1999). As these authors also included more reliable density tracers, we have adopted their value for n<sub>H</sub> and the IR derived value for G<sub>0</sub>.

In RCW 49 we will use the values given by Y. Sheffer et al. (2014, in preparation) for five out of the six IRS fields. They model the H<sub>2</sub> rotational line emission in a manner similar to that in Sheffer et al. (2011). The IRS SH observations cited in Section 3.1 are used to derive the H<sub>2</sub> 0–0 S(1) to S(4). They additionally use LL observations to cover the H<sub>2</sub> 0–0 S(0) line. PDR models are used to fit each field individually, giving H<sub>2</sub> column densities and reference values of n<sub>0</sub>. For a typical PDR, hydrogen in the surface layer is predominantly atomic and we have n<sub>H</sub> = n<sub>0</sub>. In dense PDRs, however, the H/H<sub>2</sub> transition is pulled to the surface and the atomic fraction varies from 1 to 10<sup>-3</sup>. For our sources, this seems to apply in RCW 49 regions 2 and 4. For these two positions, we hence adopted the atomic H fraction at A<sub>V</sub> = 0.5 in the detailed models of Y. Sheffer et al. (2014, in preparation); f<sub>H</sub> = 10<sup>-3</sup> and 10<sup>-1</sup> for regions 2 and 4, respectively. It should be noted that these values are very uncertain.

The value derived for NGC 7023 by Berné & Tielens (2012) is of n<sub>H</sub> = 150 ± 100 cm<sup>-3</sup>. For Orion's Veil, we estimate the density by using the electron densities derived from [S II] optical lines by Rubin et al. (2011). From their values we use a factor 20 to derive the neutral density (Tielens 2005). In the Horsehead nebula, Abergel (2003) estimate the density just behind the ionization front to be 8 × 10<sup>3</sup> cm<sup>-3</sup>. Finally, in the case of IC 63, we consider the results from Thi et al. (2009), who give a range from 1000–5000 cm<sup>-3</sup> for the PDR density.

For all the PDRs considered here, the values of n<sub>H</sub> carry uncertainties that are difficult to estimate well. Since we are using a single value as representative for each region, we have that the variations of G<sub>0</sub>/n<sub>H</sub> within each field will be representative of variations of G<sub>0</sub> only.

### 4.5. Variations in Abundance

In this section we present results of the variations in abundance of both PAHs and C<sub>60</sub>. We calculate the abundances from the ratio of intensity of PAHs and C<sub>60</sub> with respect to G<sub>0</sub>, which in turn is related to the FIR intensity of large dust grains. Assuming that PAHs, C<sub>60</sub>, and dust compete for the same photons, we can write (Tielens 2005, Section 6.7)

$$f_C = 0.23 \left( \frac{7 \times 10^{-18} \text{cm}^2}{\sigma_{\text{FUV}}} \right) \frac{f_{\text{IR}}}{1 - f_{\text{IR}}}, \quad (2)$$

where f<sub>IR</sub> is the ratio of the PAHs or C<sub>60</sub> total intensities to the total IR intensity (Berné & Tielens 2012, given by Equation (1) for the FIR and adding the contributions from PAHs and C<sub>60</sub>, following), σ<sub>FUV</sub> is the FUV absorption cross section per carbon atom of the considered species. Finally, f<sub>C</sub> is the fraction of elemental carbon locked in the species considered, which is in turn related to the total abundance.

In this calculation, we consider the total PAH intensity as the addition of all the PAH bands below 15 µm. For NGC 2023 and Ced 201 this includes features starting from λ = 5 µm.
<figure><img src="image_11.png" /><figcaption>**Figure 11.** Variations of the carbon fraction locked in PAHs (top) and $C_{60}$ (bottom) with respect to $G_0/n_H$. The different symbols indicate the different PDRs included in this work for which $C_{60}$ is detected. For RCW 49, the numbers correspond to the individual regions considered.</figcaption></figure>

In the case of RCW 49, NGC 7023, and Orion's Veil, our coverage is limited to the range between 10–15 µm since the data in these cases have been done using the SH mode. Considering the contribution of this range to the total PAH intensity in NGC 2023 and Ced 201, we find that it is fairly constant fraction at ~0.2. We will use this value as a correction factor in order to estimate the total PAH abundance in RCW 49, NGC 7023, and Orion's Veil. A final correction is needed also for the SH data arising from the difference in flux measurements with respect of the LL data. Using the flux of the 11.2 µm feature from NGC 7023 in both the LL and SH modes we find this additional factor to be ~3. Given the mismatch between the resolution of our $G_0$ maps and the $C_{60}$ and PAH intensity maps, we re-project and convolve the intensity maps to the PACS 160 µm resolution, which is the same as the resolution of our $G_0$ maps.

In order to calculate the total intensity of $C_{60}$ and PAHs at the position of the corresponding PDRs, we measure $I_{18.9}$, $I_{PAH}$, and $I_{FIR}$ at the position were the combined emission from the available PAH bands peaks. We need to include the contributions to $I_{C_{60}}$ from the bands at 17.4 µm (which also has contributions from PAHs), 8.5 µm, and 7.0 µm which we have not detected separately. Considering a first-order approximation, we use the value of 0.6 given by Peeters et al. (2012) for the ratio of the 17.4 and the 18.9 µm band intensities, and 0.4 for the ratio of both 7.0 and 8.5 µm with respect to 18.9 µm (Berné & Tielens 2012). The exact value of these ratios, however, is not precisely determined (Bernard-Salas et al. 2012), but even when considering other calculated values for the intrinsic rates, we find differences of at most a factor of two.

Using the carbon fractions derived in the previous section, we compare the results for our sources with other values found

in the literature. In NGC 7023, the observation of growth in $C_{60}$ ranges from $10^{-5}$ to $10^{-4}$ (Berné & Tielens 2012). A study of C-rich, H-poor PNe by García-Hernández et al. (2011) shows a range in the carbon fraction in $C_{60}$ from $3 \times 10^{-5}$ up to $3 \times 10^{-3}$, with most objects falling in the range of ~$10^{-4}$. All these previous estimations include the contributions to $I_{C_{60}}$ from the bands at 17.4 µm (which also has contributions from PAHs), 8.5 µm, and 7.0 µm which we have not detected separately. The value of all the regions considered in this study range fall in the $f_C = 3 \times 10^{-5}$–$6 \times 10^{-4}$, which is in good agreement with the values found in the previously mentioned studies.

The results of $f_C$ for both PAHs and $C_{60}$ in the different PDRs with respect to $G_0/n_H$ are presented in Figure 11. For PAHs, we observe a decrease with increasing $G_0/n_H$, for $G_0/n_H < 1$, from $f_C \sim 0.08$ to $f_C \sim 0.01$. For $G_0/n_H \sim 1–4$ there is substantial scatter, even within single sources. In contrast, the case of NGC 7023 shows a clear trend, with higher abundance than in other regions in the range $10 < G_0/n_H < 100$. The behavior observed in NGC 7023 appears to be the same as that observed for the other regions in the range of $G_0/n_H < 1$.

For $C_{60}$ we do not observe a general trend that can be applied to all the regions. However, within individual objects we observe that for $G_0/n_H > 1$ there is an increase in the abundance for RCW 49, NGC 7023 and Orion's Veil. In the case of RCW 49, the increase in abundance does not apply for the full range. For $G_0/n_H < 1$ we find that the abundance of $C_{60}$ within RCW 49 appears to remain constant. Furthermore, even for $G_0/n_H > 1$ the trend is not clear, particularly considering that it includes regions 2 and 4, which have poorly determined hydrogen densities. In NGC 7023 and Orion's Veil, the increase holds for the full probed range, although for none of these
<figure><img src="image_12.png" /><figcaption>**Figure 12.** Variations of the carbon fraction locked in PAHs (top) and C<sub>60</sub> (bottom) with respect to G<sub>0</sub>. The different symbols indicate the different PDRs included in this work for which C<sub>60</sub> is detected. For RCW 49, the numbers correspond to the individual regions considered.</figcaption></figure>

regions we have points with G₀/n<sub>H</sub> < 1. We cannot conclude if this behavior applies as well for NGC 2023 or Ced 201 since we have single points for these PDRs.

Figure 12 shows again $f_C$ variations of both species, in this case with respect to $G_0$. For the PAHs we observe a decrease in abundance in NGC 7023 and, to a lesser extent, in RCW 49. The abundance in Orion's Veil is consistent with a fairly constant value or a slight increase. The case of C<sub>60</sub> abundance again shows trends within individual regions. In NGC 7023 and Orion's Veil there is an increase in C<sub>60</sub> abundance along the probed range. In Orion's Veil this increase is observed at $G_0$ > $10^4$, while for NGC 7023 the increase seems to halt when $G_0$ approaches $10^4$. On the other hand, RCW 49 shows a decrease in C<sub>60</sub> abundance, which also begins at $G_0 \sim 10^4$.

## 5. DISCUSSION

Perusing Figures 11 and 12, we conclude that within NGC 7023, there is a clear trend for a decrease in PAH abundance and an increase in C<sub>60</sub> abundance with increasing G<sub>0</sub>/n<sub>H</sub>. The other regions show a more mixed behavior. For RCW 49, for example, the PAH abundance decreases with increasing G<sub>0</sub>/n<sub>H</sub> or G<sub>0</sub>, but the trend in the C<sub>60</sub> abundance is less than clear, while for Orion's Veil the opposite is true. We attribute the clearness of the trend in NGC 7023 to the fact that the PDR is seen edge-on, while for other regions the morphological and geometrical characteristics are more uncertain and might play a role in the observed trends or lack thereof. Combining all sources together there is no obvious trend whatsoever in either species or with either variable. From this, we conclude that besides the

physical conditions—G<sub>0</sub> and n<sub>H</sub>—there must be other parameter(s) influencing the evolution of the PAHs and fullerenes.

The clear trends in the PAH and fullerene abundances with G<sub>0</sub> in NGC 7023 have been interpreted in a simple model, describing the chemical evolution of PAHs under the influence of the stellar UV photons through H loss to graphene sheets resulting ultimately, on the one hand, into small hydrocarbons due to C<sub>2</sub> fragmentation and, on the other hand to fullerenes through isomerization (Berné & Tielens 2012). We note that time is an additional factor entering such chemical models. Specifically, H loss is expected to be described by a balance between collisional hydrogenation and UV-driven dehydrogenation and thus regulated by G<sub>0</sub> and n<sub>H</sub>. However, C<sub>2</sub> loss is likely to be a time-dependent process controlled by G<sub>0</sub>. Hence the absence of a general trend in the PAH and fullerene abundances across many sources can be seen as a consequence of this time-dependence. Moreover, PAH abundance appears to be related to the age of the PDRs. We observe the highest PAH abundance in Ced 201 and NGC 7023, the youngest sources in our sample, with respective ages of 1500 and 10<sup>5</sup> yr (Aleanian et al. 2008). For the remaining PDRs we find similar, overlapping age estimates: 1–7 Myr for NGC 2023 (López-García et al. 2013), 2–3 Myr for RCW 49 (Piatti et al. 1998), and 1–3 Myr for Orion’s Veil (Flaccomio et al. 2003). For these PDRs we also find lower abundance of PAHs than in Ced 201 or NGC 7023, which is consistent with a time-dependent destruction of PAHs.

Recently, the first step in the processing from PAHs to fullerenes has been modeled in detail by Montillaud et al. (2013). The results of this model are shown in Figure 13. On the left-hand side of this figure, PAHs are fully hydrogenated while
<figure><img src="image_13.png" /><figcaption>**Figure 13.** Hydrogenation state of circumcoronene ($C_{54}H_{18}$; black) and circumovalene ($C_{66}H_{20}$; red) as a function of $G_0$ and the atomic hydrogen density. The labels of the lines indicate the fraction of the respective molecules that is fully dehydrogenated. The full circles correspond to the PDRs where we detect $C_{60}$, and the symbol size is proportional to the $C_{60}$-to-PAH ratio. The sources for which we can only derive upper limits are marked by open circles. The error bars mark the ranges for the physical conditions. (Figure adapted from Montillaud et al. 2013).</figcaption></figure>

(A color version of this figure is available in the online journal.)

the right-hand-side corresponds to completely dehydrogenated PAHs—i.e., graphene flakes. More specific, the lines in this figure represent the loci at which the hydrogenation-balance leads to a constant H-fraction on a PAH of a given size. So, the further to the right, the more of the initial PAH size distribution will be destroyed but as long as the PAH size distribution extends to ~100 C atoms, some fraction of the PAHs will remain hydrogenated. The different sources in our sample are indicated by symbols in this figure. We note that the sources without detectable $C_{60}$ are to the left of the “stability-line” for PAHs in excess of 60 C atoms, leading some credence to a photochemical transformation model of PAHs into fullerene in the ISM. However, the $C_{60}$/PAH abundance ratio does not correlate well with the position into the zone where PAHs of the relevant sizes are expected to be fully dehydrogenated, as can be seen in Figure 13. This may merely indicate that time is of the essence or that there are real fluctuations in the PAH family from region to region. Ced 201 forms an interesting exception to this as it has $C_{60}$ despite being in the “stable-zone” for PAHs with more than 60 C atoms. This is, however, a very peculiar source as it is a chance encounter of the illuminating star with a cloud. This star-cloud interaction seems to have led to a shock wave (Cesarsky et al. 2000) located in the position where the $C_{60}$ abundance is particularly high. Hence, it is tempting to speculate that in this particular source $C_{60}$ is formed by shock processing of PAHs and/or dust grains. The inferred shock velocity is quite low (~10 km s⁻¹ Witt et al. 1987) and such low velocity shocks are not expected to lead to much PAH destruction (Micelotta et al. 2010). On the other hand, the observed $C_{60}$/PAH abundance is only 0.005, and hence little processing is required.

Finally, small HAC grains have also been considered as photochemical starting points for $C_{60}$ formation in the ISM. The general processes involved are the same: the first step requires dehydrogenation of the HAC particles followed by a shrinking of the resulting cage (Micelotta et al. 2012). As for PAHs, the hydrogenation of HAC is set by a balance between reactions with H and photochemical H loss (Mennella et al.

2001, 2002; Mennella 2006) and therefore controlled by $G_0/n_H$ (Chiar et al. 2013). The main difference with the PAH model resides in the timescale for C loss. Based upon molecular dynamics calculations (Zheng et al. 2007), this C loss timescale is calculated to equal the IR-cooling timescale at 100–200 K for HACs and such a low temperature is readily attained by HAC grains containing ≈500 C atoms (or less). However, the adopted C-binding energies are only 0.36 eV (Micelotta et al. 2012) and this is much less than expected for carbon cages. As emphasized in the original molecular dynamics study, these results may not be statistically significant (S. Irle, private communication; Zheng et al. 2007). Further laboratory studies will have to settle this issue.

## 6. CONCLUSIONS

We present a survey of $C_{60}$ in PDRs. While in NGC 2023 the presence of fullerenes had already been established (Sellgren et al. 2010; Peeters et al. 2012), the detection on RCW 49 and Ced 201 represent new sources where $C_{60}$ is confirmed to be present. We also quantified the abundances as fraction of elemental carbon and found values consistent with other studies of $C_{60}$ in PDRs, with values ranging from ~3 × 10⁻⁵ to ~6 × 10⁻⁴. Furthermore, the values we derive for $T_d$ from FIR observations indicate that these regions have temperatures ranging from 20 to 40 K, which is too low for $C_{60}$ to be emitting in grains, and supports the idea of a gas phase species or small sized clusters undergoing stochastic heating (Sellgren et al. 2010). The well-known strong IR bands at 6.2, 7.7, 8.6, 11.2, 12.7, and 16.4 µm, show a very similar spatial behavior, with minor variations in all sources. In contrast, the spatial distribution of the 18.9 µm band is very different and we conclude that this band has a different carrier (i.e., $C_{60}$ Cami et al. 2010; Sellgren et al. 2010) than the other bands (i.e., PAHs).

While some of the sources appear to have trends, we find no relation between $C_{60}$ or PAH abundance and either $G_0/n_H$ or $G_0$ when considering all PDRs together. We consider age as a factor explaining the lack of a general trend. Comparing our observational results to the model predictions of Montillaud et al. (2013) for PAH dehydrogenation, we find that regions where $C_{60}$ is detected have physical conditions consistent with full dehydrogenation of PAHs of at least 60 C atoms, with the exception of Ced 201. The conditions of regions where only upper limits for $C_{60}$ abundance could be derived, support only partial dehydrogenation of PAHs of the relevant size. These results support models where the dehydrogenation of carbonaceous species is the first step toward $C_{60}$ formation (Berné & Tielens 2012; Micelotta et al. 2012).

More observations aimed at measuring variations of $C_{60}$ abundance with respect to both PAHs and HAC, as well as better determinations of $n_H$ are needed to confirm the reality of these trends and the more likely parent species. Better constraints on the age estimates and the study of additional PDRs with significant age differences are needed to test our hypothesis about the effect of age in generating the seemingly independent trends for each PDR.

Studies of interstellar chemistry at Leiden Observatory are supported through advanced-ERC grant 246976 from the European Research Council, through a grant by the Dutch Science Agency, NWO, as part of the Dutch Astrochemistry Network, and through the Spinoza premie from the Dutch
Science Agency, NWO. P.C. acknowledges support from a Huygens fellowship at Leiden University.

We thank B. Ochsendorf and H. Andrews for facilitating us their reduced spectra of B 33 and IC 63, respectively. We also thank J. Montillaud for the data on the dehydrogenation state of PAHs. We also acknowledge the anonymous referee, whose comments helped us in improving this paper.

This work is based in part on observations made with the *Spitzer Space Telescope*, which is operated by the Jet Propulsion Laboratory, California Institute of Technology under a contract with NASA.

The Herschel spacecraft was designed, built, tested, and launched under a contract to ESA managed by the Herschel/Planck Project team by an industrial consortium under the overall responsibility of the prime contractor Thales Alenia Space (Cannes), and including Astrium (Friedrichshafen) responsible for the payload module and for system testing at spacecraft level, Thales Alenia Space (Turin) responsible for the service module, and Astrium (Toulouse) responsible for the telescope, with in excess of 100 subcontractors.

PACS has been developed by a consortium of institutes led by MPE (Germany) and including UVIE (Austria); KU Leuven, CSL, IMEC (Belgium); CEA, LAM (France); MPIA (Germany); INAF-IFSI/OAA/OAP/OAT, LENS, SISSA (Italy); and IAC (Spain). This development has been supported by the funding agencies BMVIT (Austria), ESA-PRODEX (Belgium), CEA/CNES (France), DLR (Germany), ASI/INAF (Italy), and CICYT/MCYT (Spain).

## REFERENCES

Abergel, A., Teyssier, D., Bernard, J. P., et al. 2003, A&A, 410, 577

Alecian, E., Catala, C., Wade, G. A., et al. 2008, MNRAS, 385, 391

Bernard-Salas, J., Cami, J., Peeters, E., et al. 2012, ApJ, 757, 41

Berné, O., Mulas, G., & Joblin, C. 2013, A&A, 550, L4

Berné, O., & Tielens, A. G. G. M. 2012, PNAS, 109, 401

Boersma, C., Rubin, R. H., & Allamandola, L. J. 2012, ApJ, 753, 168

Brand, J., & Blitz, L. 1993, A&A, 275, 67

Brown, A. G. A., de Geus, E. J., & de Zeeuw, P. T. 1994, A&A, 289, 101

Cami, J., Bernard-Salas, J., Peeters, E., & Malek, S. E. 2010, Sci, 329, 1180

Carraro, G., Turner, D., Majaess, D., & Baume, G. 2013, A&A, 555, A50

Casey, S. C. 1991, ApJ, 371, 183

Cesarsky, D., Lequeux, J., Ryter, C., & Gérin, M. 2000, A&A, 354, L87

Chase, B., Herron, N., & Holler, E. 1992, JPhCh, 96, 4262

Cherchneff, I., Le Teuff, Y. H., Williams, P. M., & Tielens, A. G. G. M. 2000, A&A, 357, 572

Chiar, J. E., Tielens, A. G. G. M., Adamson, A. J., & Ricca, A. 2013, ApJ, 770, 78

Choi, C. H., Kertesz, M., & Mihaly, L. 2000, JPCA, 104, 102

Churchwell, E., Whitney, B. A., Babler, B. L., et al. 2004, ApJS, 154, 322

Draine, B. T. 2003, ARA&A, 41, 241

Dupac, X., del Burgo, C., Bernard, J.-P., et al. 2003, A&A, 404, L11

Fabian, J. 1996, PhRvB, 53, 13864

Flaccomio, E., Damiani, F., Micela, G., et al. 2003, ApJ, 582, 398

Foing, B. H., & Ehrenfreund, P. 1994, Natur, 369, 296

Fuente, A., Martín-Pintado, J., & Gaume, R. 1995, ApJL, 442, L33

García-Hernández, D. A., Iglesias-Groth, S., Acosta-Pulido, J. A., et al. 2011, ApJL, 737, L30

Gielen, C., Cami, J., Bouwman, J., Peeters, E., & Min, M. 2011, A&A, 536, 54

Habing, H. J. 1968, BAN, 19, 421

Harvey, P. M., Thronson, H. A., Jr., & Gatley, I. 1980, ApJ, 235, 894

Hög, E., Fabricius, C., Jr., Makarov, V. V., et al. 2000, A&A, 355, L27

Houck, J. R., Roellig, T. L., van Cleve, J., et al. 2004, ApJS, 154, 18

Iglesias-Groth, S., Cataldo, F., & Manchado, A. 2011, MNRAS, 413, 213

Ijima, S. 1991, Natur, 354, 56

Jansen, R. A., van Dishoeck, E. F., & Black, J. H. 1994, A&A, 282, 605

Kemper, C., Spaans, M., Jansen, D. J., et al. 1999, ApJ, 515, 649

Knapp, S. L., Brown, R. L., & Kuiper, T. B. H. 1975, ApJ, 196, 167

Krätschmer, W., Lamb, L. D., Fostiropoulos, K., & Huffman, D. R. 1990, Natur, 347, 354

Kroto, H. W., Heath, J. R., O'Brien, S. C., Curl, R. F., & Smalley, R. E. 1985, Natur, 318, 162

Lada, E. A., DePoy, D. L., Evans, N. J., II, & Gatley, I. 1991, ApJ, 371, 171

López-García, M. A., López-Santiago, J., Albacete-Colombo, J. F., Pérez-González, P. G., & de Castro, E. 2013, MNRAS, 429, 775

Maier, J. P. 1994, Natur, 370, 423

Menéndez, J., & Page, J. B. 2000, in Light Scattering in Solids VIII: Fullerenes, Semiconductor Surfaces, Coherent Phonons, ed. M. Cardona, Güntherodt (Berlin: Springer), 27

Mennella, V. 2006, ApJL, 647, L49

Mennella, V., Brucato, J. R., Colangeli, L., & Palumbo, P. 2002, ApJ, 569, 531

Mennella, V., Muñoz Caro, G. M., Ruiterkamp, R., et al. 2001, A&A, 367, 355

Micelotta, E. R., Jones, A. P., Cami, J., et al. 2012, ApJ, 761, 35

Micelotta, E. R., Jones, A. P., & Tielens, A. G. G. M. 2010, A&A, 510, A36

Moffat, A. F. J., Shara, M. M., & Potter, M. 1991, AJ, 102, 642

Montillaud, J., Joblin, C., & Toublanc, D. 2013, A&A, 552, 15

Mookerjea, B., Sandell, G., Jarrett, T. H., & McMullin, J. P. 2009, A&A, 507, 1485

Ohama, A., Dawson, J. R., Furukawa, N., et al. 2010, ApJ, 709, 975

Otsuka, M., Kemper, F., Sargent, B., et al. 2012, in ASP Conf. Ser. 458, Galactic Archaeology: Near-Field Cosmology and the Formation of the Milky Way, ed. W. Aoki, M. Ishigaki, T. Suda, T. Tsujimoto, & N. Arimoto (San Francisco, CA: ASP), 137

Pankonin, V., & Walmsley, C. M. 1976, A&A, 48, 314

Peeters, E., Tielens, A. G. G. M., Allamandola, L. J., & Wolfire, M. G. 2012, ApJ, 747, 44

Piatti, A. E., Bica, E., & Claria, J. J. 1998, A&AS, 127, 423

Pilbratt, G. L., Riedinger, J. R., Passvogel, T., et al. 2010, A&A, 518, L1

Pilleri, P., Montillaud, J., Berné, O., & Joblin, C. 2012, A&A, 542, A69

Poglitsch, A., Waelkens, C., Geis, N., et al. 2010, A&A, 518, L2

Rauw, G., De Becker, M., Nazé, Y., et al. 2004, A&A, 420, L9

Rauw, G., Manfroid, J., Gosset, E., et al. 2007, A&A, 463, 981

Roberts, K. R. G., Smith, K. T., & Sarre, P. J. 2012, MNRAS, 421, 3277

Rubin, R. H., Simpson, J. P., O'Dell, C. R., et al. 2011, MNRAS, 410, 1320

Sellgren, K., Werner, M. W., Ingalls, J. G., et al. 2010, ApJL, 722, L54

Sheffer, Y., Wolfire, M. G., Hollenbach, D. J., Kaufman, M. J., & Cordier, M. 2011, ApJ, 741, 45

Shetty, R., Kauffmann, J., Schnee, S., & Goodman, A. A. 2009, ApJ, 696, 676

Smith, J. D. T., Armus, L., Dale, D. A., et al. 2007a, PASP, 119, 1133

Smith, J. D. T., Draine, B. T., Dale, D. A., et al. 2007b, ApJ, 656, 770

Steiman-Cameron, T. Y., Haas, M. R., Tielens, A. G. G. M., & Burton, M. G. 1997, ApJ, 478, 261

Thi, W.-F., van Dishoeck, E. F., Bell, T., Viti, S., & Black, J. 2009, MNRAS, 400, 662

Tielens, A. G. G. M. 2005, The Physics and Chemistry of the Interstellar Medium (1st ed.; Cambridge: Cambridge Univ. Press)

Tielens, A. G. G. M. 2008, ARA&A, 46, 289

Werner, M. W., Roellig, T. L., Low, F. J., et al. 2004, ApJS, 154, 1

Westerlund, B. 1960, ArA, 2, 419

Whiteoak, J. B. Z., & Uchida, K. I. 1997, A&A, 317, 563

Witt, A. N., Graff, S. M., Bohlin, R. C., & Stecher, T. P. 1987, ApJ, 321, 912

Young Owl, R. C., Meixner, M. M., et al. 2002, ApJ, 578, 885

Zhang, Y., & Kwok, S. 2011, ApJ, 730, 126

Zheng, G., Wang, Z., Irle, S., & Morokuma, K. 2007, J. Nanosci. Nanotechnol., 7, 1662