text:
# ISOCAM-CVF 5-12𝜇m spectroscopy of Ultraluminous Infrared Galaxies 1footnote 11footnote 1Based on observations with ISO, an ESA project with instruments funded by ESA Member States (especially the PI countries: France, Germany, the Netherlands and the United Kingdom) with the participation of ISAS and NASA

## 1 Introduction

Recently, mid infrared spectroscopy has been used increasingly as a tool for studies of ultraluminous infrared galaxies (ULIRGs, $`L_{\mathrm{IR}}>10^{12}L_{\mathrm{}}`$, see Sanders & Mirabel (1996) for a review). The presence of starbursts and active galactic nuclei (AGN) in ULIRGs (and coexistence in some of them) has been known for quite some time, but the question which of the two dominates the luminosity has been difficult to answer because of the large columns of obscuring dust found towards the nuclear regions of these gas- and dust-rich systems. High sensitivity and complete coverage of the infrared spectrum by ESA’s Infrared Space Observatory (ISO, Kessler et al. 1996) have considerably advanced the use of infrared spectroscopy to penetrate this obscuring dust, and to probe for the sources of the huge luminosity of ULIRGs.

Fine structure line and aromatic emission feature observations with ISO-SWS and ISOPHOT-S of a sample of 15 bright ULIRGs suggest that most ULIRGs are predominantly starburst-powered (Genzel et al., 1998). The PAH method has been extended to a larger sample (Lutz et al., 1998; Rigopoulou et al., 1999), allowing to probe for evolutionary effects expected in a scenario where starburst activity gives way to a quasar-like AGN buried inside the ULIRG (Sanders et al., 1988). Comparison of optical and mid-infrared spectroscopic diagnostics (Lutz et al., 1999; Taniguchi et al., 1999) demonstrated surprisingly good agreement if optical LINERs are interpreted as starbursts. This finding also suggests that AGN in ULIRGs usually make their presence known optically at least in certain directions, instead of being fully embedded by large obscuring columns of dust.

This paper presents the result of a solicited program in ISO open time, which addresses two main issues only partly covered by previous work. First, and most important, at which luminosity does the transition from ‘predominantly starburst-powered’ to ‘predominantly AGN-powered’ occur? On one hand, the ISO observations of Genzel et al. (1998) and Lutz et al. (1998) suggest that not only luminous infrared galaxies with $`L_{\mathrm{IR}}<10^{12}L_{\mathrm{}}`$, but also most ULIRGs above this threshold are predominantly starburst powered. On the other hand, optical and mid-infrared observations of hyperluminous infrared galaxies (Hines et al., 1995, 1999; Taniguchi et al., 1997; Aussel et al., 1998) support a previous consensus that AGN dominate these hyperluminous systems (but note the intriguing possibility that far-infrared peaks in some luminous AGN may trace even higher star formation related luminosities, e.g. Haas et al. 1998). The luminosity at which local ULIRGs ($`z<0.4`$) – on average – switch from starburst to AGN dominated is a quantity with implications for various fields ranging from the interpretation of sources found in recent submm surveys (Hughes et al., 1998; Barger et al., 1998) to explanation of the hard X-ray background (Comastri et al., 1995). Simple arguments based on the free fall time of the gas mass concentrated in the inner region of ULIRGs and the nucleosynthesis efficiency suggest a maximum luminosity of a starburst event approaching $`10^{13}L_{\mathrm{}}`$ (Heckman, 1994), the limit depending on gas mass, spatial scale, and initial mass function properties. Observational evidence is clearly needed. Despite its size of $``$60 sources, the ULIRG sample presented by Lutz et al. (1998) and Rigopoulou et al. (1999) is not optimal for this task: It is selected basically as flux limited at 60$`\mu `$m and dominated by low to moderate luminosity ULIRGs, because of the steep ULIRG luminosity function (Sanders & Mirabel, 1996). The preliminary analysis of Lutz et al. (1998) in the mid-infrared and a recent study of Veilleux et al. (1999) in the optical/NIR infrared showed that the presence of AGN in ULIRGs increases with the luminosity of the object. However, these studies clearly call for improved statistics above $`10^{12.3}L_{\mathrm{}}`$ in order to better determine the transition to AGN-dominated systems.

The second issue addressed by this paper is related to the use of PAH features as an AGN/starburst diagnostic. Mid-infrared spectra of most galaxies show the 6.2, 7.7, 8.6, and 11.3$`\mu `$m features attributed to aromatic carbonaceous material (Duley & Williams, 1981; Léger & Puget, 1984; Sakata et al., 1987; Papoular et al., 1989). Among several popular designations for these bands will be designated hereafter by the terms PAH (polycylic aromatic hydrocarbon) or UIB (unidentified infrared bands, because of remaining uncertainties on the precise nature of the aromatic carrier). Groundbased observations of these features and a companion at 3.3$`\mu `$m first demonstrated that their equivalent width is larger in starburst galaxies than in classical AGNs (Moorwood, 1986; Roche et al., 1991). ISO spectroscopy has further strengthened this link by demonstrating the anti-correlation between feature strength relative to the continuum and the ionization state of the gas (Genzel et al., 1998). Spatially resolved ISOCAM-CVF observations of nearby AGN support this interpretation by showing the PAH features contrast to be weak near the central AGN but strong in the circumnuclear region likely dominated by star formation (e.g. Cen A, NGC 1068: Laurent et al. (1999), Circinus: Moorwood (1999)). Despite this solid empirical basis, special care is required in analyzing PAH spectra of ULIRGs for AGN or starburst dominance. All components of the spectrum – continuum and emission features, AGN or starburst – may be unusually obscured, and the limited wavelength coverage of ISOPHOT-SL spectra (5.8-11.6$`\mu `$m observed wavelength) makes continuum definition difficult. Here, the extended wavelength range of ISOCAM-CVF (5-16.5$`\mu `$m) is beneficial in reaching the long wavelength side of the broad 9.7$`\mu `$m silicate absorption feature.

Our paper is organized as follows: We discuss sample selection, observations and data analysis in Section 2. Section 3 presents the observational results. Section 4 establishes a new quantitative method to determine the contribution of starburst and AGN activity to the mid-infrared spectra, and compares it with the method used by Genzel et al. (1998), Lutz et al. (1998) and Rigopoulou et al. (1999). In Section 5, we discuss the properties of ULIRGs as a function of luminosity, combining our sample with the results of the ISOPHOT-S sample and with ISOCAM results on HYLIRGs. Finally, we conclude in section 6.

## 2 Observations and data reduction

### 2.1 Sample and observing strategy

The solicited observing program presented here was defined in the last third of the ISO mission when it had become apparent that low resolution spectroscopy of ULIRGs was a promising tool, and that studies of ULIRGs needed to be extended towards higher luminosities. The sample selection is strongly driven by the visibility constraints of the ISO satellite during the final six months of the mission. ULIRGs were selected from the samples of Fisher et al. (1995), Kim (1995), Clements et al. (1996), and the QDOT sample (Lawrence et al. 1999). Six sources were chosen from these catalogues with $`L_{\mathrm{IR}}>10^{12.5}L_{\mathrm{}}`$ that had reasonable ISO visibility during the remaining mission, plus an additional ten slightly lower luminosity sources that could be accommodated within the available observing time. Table 1 lists basic properties of the sample galaxies. Throughout this paper, we follow the convention reported by Sanders & Mirabel (1996) to compute the 8-1000$`\mu `$m luminosity $`L_{\mathrm{IR}}`$, adopting H0=75 km s-1 Mpc-1 and q0=0.5, and IRAS FSC fluxes ( IRAS 03521+0028, IRAS 18030+0705: PSC). Most of our sources are lacking detections in some of the IRAS photometric bands. In those cases, we replaced the missing fluxes by estimates based on the 60$`\mu `$m fluxes and average colors of the Sanders et al. (1988) ULIRGs. The observations included the source IRAS 03452+2320 taken from the sample of Fisher et al. (1995), which is not listed in Table 1. We did not detect a signal at its position. We believe this reflects an infrared cirrus related misidentification in the Fisher et al. (1995) sample rather than mid-infrared faintness of a true ULIRG, since Crawford et al. (1996) failed to detect the same source in radio continuum.

Our sample is certainly heterogeneous : it was selected from ULIRG surveys of different limiting 60$`\mu `$m flux, with a preference for high luminosity targets but not implementing a strict lower cutoff, and strongly driven by ISO visibilities. None of the selection criteria, however, should significantly bias our sample in mid-IR spectral properties or AGN content. In particular, no AGN-related IRAS color criteria like the 25/60$`\mu `$m flux ratio were applied (de Grijp et al., 1985; Sanders et al. 1988a, ). A slight tendency towards favoring AGN at higher redshifts/luminosities is induced by the 60$`\mu `$m selection which is effectively at $``$45$`\mu `$m rest wavelength for the most distant targets. We have not attempted to correct for this since detailed far-infrared spectral energy distributions are unavailable for these faint IRAS sources, and since the effect will be negligible for our small and still fairly local sample. The sample is hence suited to study the AGN content of ULIRGs, bridging the gap between the lower luminosity ULIRGs studied with ISOPHOT-S (Lutz et al., 1998; Rigopoulou et al., 1999) and the two HyperLuminous Infrared Galaxies (HYLIRGs) observed by ISO Taniguchi et al. (1997); Aussel et al. (1998). Figure ISOCAM-CVF 5-12$`\mu `$m spectroscopy of Ultraluminous Infrared Galaxies 111Based on observations with ISO, an ESA project with instruments funded by ESA Member States (especially the PI countries: France, Germany, the Netherlands and the United Kingdom) with the participation of ISAS and NASA shows the distribution of infrared luminosities both for our sample and the combined ISO database.

The observations presented here used the Circular Variable Filter (CVF) mode of the ISOCAM camera (Cesarsky et al., 1996) rather than the ISOPHOT-S mode with which most other low resolution ULIRG spectra were obtained and which is efficient in building a large sample of relatively bright sources. The ultimate sensitivity achieved in a long CVF integration time is better, and the need for stabilization overheads was no major concern for our moderate number of very faint sources. In addition, the extended long wavelength coverage is crucial in targeting objects up to z$``$0.4, and in probing an extended rest wavelength range. Since our sources are expected to be spatially unresolved, we have chosen the most sensitive setting with 6″ pixels and a total FOV of about 3′. The observed wavelength range varies from source to source, the short end being defined by 5$`\mu `$m rest wavelength and the long end either by 12$`\mu `$m rest wavelength or, for high redshift sources, by 16$`\mu `$m observed wavelength where the sensitivity of CVF observations starts to decrease quickly. In order to increase redundancy and to enable correction of the transients of the ISOCAM detector, each CVF scan has been done in up and down direction. To reconcile the need for detector stabilization (10 frames of 5sec exposure at each wavelength) with the total number of wavelength steps and practical limits on observing time, we have executed up- and down-scans with interleaved wavelength steps. The total observing time per source including overheads was typically 2.5 hours. Each final spectrum is composed by about 120 wavelength steps and has an average resolution of $`R35`$.

### 2.2 Data reduction

Reductions have been carried out with the Cam Interactive Analysis (CIA) software222CIA is a joint development by ESA Astrophysics division and the ISOCAM consortium led by the ISOCAM PI, C. Cesarsky, Direction des Sciences de la Matière, CEA, France.. The dark current has been estimated with the dark model described by Biviano & Sauvage (1999). This gives better results than using the two offset observations we have taken before and after the CVF observation. We used the multi-median resolution deglitching method (Starck et al., 1999) to remove glitches, and then additionally deglitched manually a square of 3 by 3 pixels centered on the source. To correct the transient effect of the detector the inversion method developed at IAS (Coulais & Abergel, 1999) was used. The calibration of the flux was performed using the standard spectral response function of CIA.

Appendix A presents a reduction scheme which allows to quantify the spectrum of a target, with its S/N ratio and the systematic error introduced by the reduction method (the total uncertainty on the flux is about 5 mJy). It also allows us to determine the position of the fairly faint source with sub-pixel accuracy. The typical uncertainty of the position is 5″, including S/N and ISO pointing uncertainty.

One target (IRAS 00406-3127) of the set had not been acquired by ISOCAM because of a telemetry drop. For IRAS F02115+0226 there was no source detection. The S/N ratio observed for IRAS F02455-2220 is not high enough to reliably determine the PSF correction described in the Appendix. However, the observed spectrum can be taken as a lower limit, using only the signal from the brightest pixel. For the remaining 13 targets, the data reduction has been applied successfully.

## 3 Results

The rest-frame spectra for the 13 sources are shown in Figure ISOCAM-CVF 5-12$`\mu `$m spectroscopy of Ultraluminous Infrared Galaxies 111Based on observations with ISO, an ESA project with instruments funded by ESA Member States (especially the PI countries: France, Germany, the Netherlands and the United Kingdom) with the participation of ISAS and NASA. The shaded confidence area around the spectrum represents the RMS plus the systematic error. Table 1 lists the position of the sources as obtained by the PSF fitting method (see Appendix), and the IRAS fluxes. All targets were found to be point like sources at ISOCAM resolution. One target has been found to be significantly offset from the IRAS PSC position: IRAS 18030+0705. The offset is about one arcminute. We checked the original IRAS data using XSCANPI and found a source with fluxes similar to those listed in the PSC, but centered on the ISOCAM position. The origin of this position discrepancy is unknown. The FSC fluxes were adopted for the ULIRGs detected by ISOCAM except for IRAS 18030+0705 and IRAS 03521+0028 for which we used PSC fluxes.

¿From previous studies, the following components are expected to be present in ULIRG spectra, but not necessarily together in one spectrum :

* UIB emission features at 6.2, 7.7, 8.6 and 11.3 microns, mainly originating from star forming regions of the galaxies. This UIB emission can be also coupled with a continuum of hot dust present in the most energetic H II regions of the starburst.
* An AGN continuum from the hot dust present in the vicinity of the AGN, e.g. in a dust torus or at larger scale in the narrow line region.
* Strong absorption, both continuous and in features like the silicate 9.7 $`\mu `$m one. Absorption may arise either in the surroundings of the AGN or the larger scale ISM. The amount of absorption may be different for AGN and for starburst.

The silicate feature at 9.7$`\mu `$m that is usually found in absorption may also mimic emission features. In dust configurations as in some torus models (Pier & Krolik, 1992), self-absorbed silicate emission can create an emission peak near 8$`\mu `$m that is superficially similar to a faint UIB emission and may be observed in one HYLIRG (Taniguchi et al., 1997). Discrimination of such a ‘fake UIB’ peak and real UIB emission has to rely on the detailed shape of the 8$`\mu `$m peak and presence of the other UIB features.

Table 2 summarizes the presence of the UIB bands. About half of the targets exhibit a spectrum with at least the 6.2, 7.7 and 11.3 $`\mu `$m mid infrared bands clearly visible. The remaining targets do not exhibits all these bands, and two none of them significantly. It is striking that almost all the observed ULIRGs show a strong feature at or near 7.7 $`\mu `$m. IRAS F23529-2919 is an exception. Here, the continuum might be so strong that it can dominate the emission of the 7.7 $`\mu `$m emission band. Also, the shape of the 7.7$`\mu `$m feature is very different in two other galaxies from the typical UIB band. IRAS 22192-3211 and IRAS 00275-2859 have a band which is broader than the one observed in a “normal” starburst galaxy. The signal to noise for IRAS 23515-2917 is lower so it is more difficult to draw a definitive conclusion. A strong continuum has been detected in the two former targets, which supports the presence of an AGN in both sources. The change of the shape of the 7.7 $`\mu `$m feature can be explained in two ways. The carriers of the UIB might be processed in the AGN environment, changing the shape of the spectral emission of these carriers. Such changes in the shape of the 7.7/8.6$`\mu `$m complex have in fact been observed in some ultra-compact H II regions (Cesarsky et al. 1996a, ). Alternatively, a “fake” 7.7 $`\mu `$m feature could arise when one observes a high optical depth source with self-absorbed silicate emission. There is no agreement in detail, however, between the broadish 8$`\mu `$m features in these ULIRGs and the spectral shapes predicted by torus models (Pier & Krolik, 1992). A further problem of a ‘torus’ identification of these 8$`\mu `$m features is that they are not observed in the spectra of template Seyferts and QSOs (Rigopoulou et al., 1999; Clavel et al., 2000).

Another source, IRAS F00183-7111, has a maximum near 8$`\mu `$m which we do not ascribe to UIBs. Its spectrum is remarkably similar to some deeply obscured objects (d’Hendecourt et al. 1996, NGC 7538 IRS1 : D.Cesarsky unpublished), showing a strong 5-8$`\mu `$m continuum and an extremely deep silicate absorption. Even an apparent “6.5 $`\mu `$m emission” is indicated in this object which is very likely the result of the absorption of water ice at 6 $`\mu `$m and an absorption feature at 6.85$`\mu `$m which is of unclear origin (Schutte, 1998). The spectrum of IRAS F00183-7111 is thus dominated by a heavily absorbed (AGN) continuum with a possible weak 7.7$`\mu `$m emission, which might be either due to a classical UIB or self-absorbed silicate emission. In total, it appears to be the only unambiguous example in our data for a type of source considered typical in the classical Sanders et al. (1988) scenario: an extremely obscured AGN.

## 4 Quantitative diagnostics of the Starburst and AGN activity in ULIRGs

The use of the line-to-continuum ratio of the 7.7 $`\mu `$m band (hereafter L/C) as proposed by Genzel et al. (1998) and Lutz et al. (1998) to distinguish in the mid-infrared between starburst-dominated and AGN-dominated ULIRGs is based on evidence that a starburst galaxy presents strong UIB features while AGN presents strong continuum with a low contrast emission features in the nuclear part of the galaxy. This is suggested by ground based observations Moorwood (1986); Roche et al. (1991), ISO spectroscopy Genzel et al. (1998); Rigopoulou et al. (1999), and ISO spectro-imaging Moorwood (1999); Laurent et al. (1999). This also implies that the UIB luminosity of a galaxy scales with its star formation. Recent studies tend to reinforce this idea : in spiral galaxies, the global H$`\alpha `$ emission of the disk is correlated with the total MIR emission of UIB Vigroux et al. (1999); Roussel et al. (2000).

When using these features as tracers of the sources of the total ULIRG luminosity, a further assumption is that the mid-infrared actually probes the active regions and not just a surface layer. This is not trivial since the active regions of ULIRGs are small Condon et al. (1991); Soifer et al. (2000) and dusty. Depending on dust configuration, a major part of the luminosity could be hidden from mid-infrared view. Genzel et al. (1998) have addressed this problem using the starburst-type mid-infrared nebular spectra of ULIRGs. For ULIRGs, the ratio of the inferred starburst ionizing luminosity to the bolometric luminosity approaches that for template starbursts, indicating that more than half of the luminosity has been traced by the mid-infrared observations. Similar analyses can be made using the luminosity of the PAH features. Rigopoulou et al. (1999) and Genzel & Cesarsky (2000) find that observed PAH luminosities of ULIRGs are well above those of lower luminosity sources, the ratio to bolometric luminosity being only a factor 1.4-2 lower than in starbursts. This factor is plausibly explained by the higher extinction in ULIRGs and perhaps the effect of the more intense radiation field on the PAH emission. We have tentatively repeated this comparison using the extinction correction derived from our fits described below, which however give only fairly uncertain extinction for individual noisy spectra. The mean discrepancy narrows down to a factor 1.2. Thus, the mid-infrared diagnostics are probing the major part of the luminosity.

For a narrow line and emission from a single region, both line and continuum are obscured by the same amount. The width of the PAH feature and the strong wavelength dependence of extinction at the onset of the silicate feature make the PAH L/C less ideal. In some objects the absorption might be so strong that it can mimic a 7.7 band. There are, however, differences in shape and the L/C will normally stay below the cutoff of 1 adopted in previous papers. IRAS F00183-7111 is the only clear example in our sample. Another technical point that has to be verified is the determination of the continuum under the 7.7$`\mu `$m feature. In the L/C diagnostic, this continuum is deduced from the linear interpolation of the continuum levels at 5.9 and 10.9$`\mu `$m, i.e. the L/C ratio is affected by the flux uncertainties at these wavelengths. This is also a systematic source of uncertainty since the 10.9$`\mu `$m point is still within the silicate feature and does not probe a clean continuum.

We have explored alternative methods in order to verify the robustness of the L/C ratio with respect to the systematic error induced in the 7.7$`\mu `$m continuum determination. First, we constructed the “narrow band” ratio between the flux density at 7.7 and 6$`\mu `$m. Strong UIB emission will create a high ratio which is not reached by a pure continuum unless it is unrealistically steep over this small wavelength interval. This diagnostic should hence identify strong UIB emission easily but will face problems discriminating weak UIBs from continuum slope variations. The strength of this diagnostic is to avoid any assumption on the shape of the spectrum in the region obscured by the silicate feature. Figure ISOCAM-CVF 5-12$`\mu `$m spectroscopy of Ultraluminous Infrared Galaxies 111Based on observations with ISO, an ESA project with instruments funded by ESA Member States (especially the PI countries: France, Germany, the Netherlands and the United Kingdom) with the participation of ISAS and NASA shows the correlation between the L/C and F(7.7)/F(6) diagnostics seen for ULIRGs and template spectra observed by PHOT-S (Rigopoulou et al., 1999), SWS (Sturm et al., 2000), and CAM-CVF (this paper). The dispersion in this plot is explained by the variation in continuum slope and the systematic error of the L/C discussed before. The good correlation of the two diagnostics supports the use of the L/C diagnostic to identify systems that are clearly starburst-like and show strong UIBs. By definition of the F7.7/F6 diagnostic, Figure ISOCAM-CVF 5-12$`\mu `$m spectroscopy of Ultraluminous Infrared Galaxies 111Based on observations with ISO, an ESA project with instruments funded by ESA Member States (especially the PI countries: France, Germany, the Netherlands and the United Kingdom) with the participation of ISAS and NASA cannot test the effect of continuum placement uncertainty in the L/C method on measurements of faint residual UIBs in AGN-like systems. The line F(7.7$`\mu `$m)/F(6$`\mu `$m)=L/C+1 is also overplotted in Figure ISOCAM-CVF 5-12$`\mu `$m spectroscopy of Ultraluminous Infrared Galaxies 111Based on observations with ISO, an ESA project with instruments funded by ESA Member States (especially the PI countries: France, Germany, the Netherlands and the United Kingdom) with the participation of ISAS and NASA to indicate the asymptotic case where the continuum is flat.

To investigate further how the different spectral components influence the MIR spectrum of a galaxy, we created another diagnostic whose basis is a more physical parameterization of the spectrum.

### 4.1 Presentation of the model

Our aim is an empirical fitting procedure using the complete mid-infrared low resolution spectra to derive the contribution of the AGN and starburst components to the mid-infrared spectrum. In a second step, results will be extrapolated to the bolometric luminosities. The fit is based on template spectra of dusty starbursts and AGN, and applied to the whole set of ULIRG and template spectra observed by PHOT-S Rigopoulou et al. (1999), SWS Sturm et al. (2000), and CAM-CVF (this paper). Since absorption plays an important role for the spectra of dust-rich ULIRGs, we incorporate the effects of absorption that is higher than in the templates. In principle, explaining the mid-infrared spectra of ULIRGs is a problem of three dimensional radiative transfer requiring a consistent treatment of dust, ices, and PAH features. Considerable insight can be gained, however, from a simple superposition of a screen-obscured starburst and a screen-obscured AGN. Arguments in favour of this approach are: (1) The mid-infrared diagnostic reaches most of the active region. Unknown more obscured sources will be minor contributors. (2) The use of an observed AGN template implicitly accounts for warm dust close to the AGN. Screen extinction is a reasonable first approximation for the effect of even more dust around the AGN, since this dust will mostly just absorb in the mid-infrared consider here, and reemit in the far-infrared. (3) An aggregate of star forming regions is in principle better represented by a case where emission and absorption are well mixed. Mixed case extinctions are difficult to determine, however, unless there are probes measured at high S/N over a large range of optical depths, a requirement not met by our spectra. The screen case adopted here will give similar results for optical depths around one, fall short at higher optical depth, but be conservative in the sense of not inferring excessive starburst activity in uncertain cases.

#### 4.1.1 The Starburst component

Genzel et al. (1998), Rigopoulou et al. (1999) and Sturm et al. (2000) showed that in the range of the 6-9$`\mu `$m UIB bands there is little variation in starburst MIR spectra observed by SWS and PHOT-S. Starburst spectra present some variations at longer wavelength, however. Vigroux et al. (1996), noted that in some regions of NGC 4038/4039 the flux at 15 $`\mu `$m can rise dramatically and dominate the mid infrared spectrum. It can be shown that the MIR spectrum of starburst galaxies is the superposition of two components : one from photodissociation regions (PDR) and one from H II regions (Tran, 1998). The PDR component is responsible for the UIB emission and the H II region for the hot dust continuum which becomes important at 12 $`\mu `$m (Cesarsky et al. 1996a, ; Verstraete et al., 1996; Roelfsema et al., 1996). This hot dust continuum is probably due to very small grains heated to temperature as high as 200K in some regions. This continuum is very faint at 7$`\mu `$m (Cesarsky et al. 1996a, ; Verstraete et al., 1996; Tran, 1998), and is not really noticeable in the spectrum of the majority of starburst galaxies where the PDR contribution (i.e. the UIB contribution) dominates in this part of the spectrum (Sturm et al., 2000). For some dwarf metal-deficient starburst galaxies like NGC 5253 and IIZw40 (Rigopoulou et al. 1999, S.Madden private communication), UIBs are nearly absent. In such a case a strong featureless continuum will not be an AGN signature. Since ULIRG are dust rich and not known to be metal-deficient this effect is likely not relevant, and we still consider a strong ULIRG continuum with weak UIBs an AGN signature.

As the (rest)-wavelength range investigated in this paper is between 5 and 12 $`\mu `$m, the contribution of the hot dust continuum from H II regions will not change drastically the spectrum of the starburst contribution. But one must keep in mind that strictly speaking a significant continuum at short wavelengths (5-8 $`\mu `$m) is not necessarily a tracer of an AGN, but could reveal a higher weight on H II regions in the ISM with respect to photodissociation regions. To break the degeneracy, one has to look at the continuum at 15 $`\mu `$m since the continuum of an H II region is steeper than the continuum of an AGN (Laurent et al., 1999).

In the present model we assume that the template spectrum of a starburst does not change in a significant way between 5 and 12 $`\mu `$m (see section 4.2 for the self consistency of this hypothesis). We have used the M82 SWS spectrum as starburst template spectrum, supposed to represent the standard line to continuum ratio of the UIB seen in a dusty starburst.

#### 4.1.2 The continuum component

Nuclear AGN spectra have low contrast UIB features and mainly present a continuum emitted by hot dust. Obviously, when looking at the total spectrum of a galaxy hosting an AGN, one naturally observes UIB bands from circumnuclear star formation or from the disk of the AGN host (see for instance (Mirabel et al., 1999; Clavel et al., 2000; Moorwood, 1999; Alexander et al., 1999)). The AGN spectrum itself can be described by an intrinsic continuum which can be well fitted by a power law333The ‘power law’ does not mean non-thermal emission, and thus it is a parameterization of the dust continuum. of variable index (Clavel et al., 2000) subject to absorption by dust especially in the region of the 9.7$`\mu `$m silicate absorption feature.

As discussed in the previous subsection, a continuum seen in the mid-infrared may also be due to H II regions. In the wavelength range considered (5-12 $`\mu `$m), there is a degeneracy between the AGN continuum and a continuum coming from an H II region, although an H II continuum is usually steeper and fainter. Our fit invokes a single ‘continuum component’ including both possibilities of an AGN continuum and an unusually strong H II continuum. We will come back to the degeneracy between these two types of continua when discussing fit results.

#### 4.1.3 Dust obscuration

Because ULIRGs are gas and dust rich galaxies (Solomon et al., 1997), the extinction is expected to play a crucial role on shaping the spectrum of these galaxies, even in the mid-infrared (see for instance Arp220, Sturm et al. (1996)). It is important to test how dust obscuration can affect both UIB dominated and featureless spectra (Rigopoulou et al., 1999). One can argue that mid-infrared ULIRG spectra could be completely dominated by absorption, and that the 7.7$`\mu `$m feature, which we identify as an UIB, could result from a deep absorption of a flat continuum (Dudley & Wynn-Williams, 1997). This could be the case in the strongly absorbed galaxy F00183-7111. We want to test more generally whether it is possible to produce an emission feature around 8 $`\mu `$m by a power law continuum affected by an absorption law or whether a simple UIB spectrum affected by absorption is preferred.

We used two absorption laws, the ‘classical’ Draine and Lee (1984) absorption law and one composite absorption law derived from the Galactic Center-SgrA\* spectrum (hereafter GC or GC-SgrA\*). The GC extinction law has a higher extinction level between 3 and 7$`\mu `$m (Lutz et al., 1996; Lutz 1999a, ). Since this law is derived from hydrogen recombination lines and does not sample the silicate feature well, we have supplemented it with points deduced from the GC continuum spectrum itself. They have been computed assuming that the intrinsic GC spectrum is a power law spectrum absorbed by a dust screen. The extinction law and the exponent of the power law have been adjusted in order to fit the extinction deduced from the hydrogen lines and the GC spectrum itself. This procedure is adequate since we cover only a limited wavelength range, and since the dust in front of the Galactic center is close to the ideal screen case. This absorption law has two main peculiarities: a flatter absorption between 3 and 7$`\mu `$m, and a narrower silicate feature than what the Draine and Lee curve predicts (see Figure ISOCAM-CVF 5-12$`\mu `$m spectroscopy of Ultraluminous Infrared Galaxies 111Based on observations with ISO, an ESA project with instruments funded by ESA Member States (especially the PI countries: France, Germany, the Netherlands and the United Kingdom) with the participation of ISAS and NASA).

These two absorption laws have been used to represent variations in dust properties as possible in extreme environments such as close to an AGN. Comparison of the M 82 and NGC 1068 SWS spectra (Sturm et al., 2000) suggests such variations at least for absorption features. Because of the possibility of such variations, absolute AV values derived from the model are somewhat uncertain. The main emphasis is to test whether the shape of the ULIRG spectra can be reproduced well.

#### 4.1.4 The different models

The two spectral components and obscuration have been used to construct four types of models in two classes:

* The first class is based on the assumption that the starburst component and the continuum component are obscured by the same dust screen. Models of this class vary in extinction properties: the dust extinction law is either a “pure” Draine and Lee law, or a “pure” GC extinction law or a mixture of the two.
* The second class allows different obscuration of the two components and is clearly more realistic. The two different contributions (continuum and starburst) are not necessarily located in the same physical region and therefore do not suffer the same extinction.

Table 3 provides the formal description of the different models used. Models 1a, 1b and 2 refer to the first class, while models 3a-d and 4 refer to the second class. In two cases (2 and 4), the extinction law is a linear combination of the GC extinction law and the Draine and Lee one. This combination is not constrained by observations but included in our modeling as a comparison to other models with different degrees of freedom.

### 4.2 Testing the validity of the model

To validate the model, the fit quality is quantified first. Second, the results from test fits to template starburst and AGN spectra are investigated, and compared to observations of “mixed” sources where the AGN and starburst relative contributions are well known from spatially resolved data.

#### 4.2.1 Quality of model fits

Figure ISOCAM-CVF 5-12$`\mu `$m spectroscopy of Ultraluminous Infrared Galaxies 111Based on observations with ISO, an ESA project with instruments funded by ESA Member States (especially the PI countries: France, Germany, the Netherlands and the United Kingdom) with the participation of ISAS and NASA shows some examples of fits to ULIRG spectra using model 3a. The figure includes AGN and starburst-like spectra of both high and low S/N. Fits using other models tend to be slightly worse but still acceptable. The first conclusion is that the models presented in Table 3 can give a reasonable fit for all ULIRG spectra. This implies that the assumptions of these models (e.g. a MIR spectrum of a ULIRG is a superposition of a continuum component and a starburst component extincted by a dust-rich ISM) are not in obvious contradiction with observations.

To better quantify the quality of the various models we adopt the indicator $`\frac{\chi ^2}{N_{\mathrm{free}}}=\frac{{\scriptscriptstyle (F_\lambda ^{\mathrm{obs}}F_\lambda ^{\mathrm{mod}})^2}}{{\scriptscriptstyle (\sigma ^{\mathrm{obs}})^2}}\frac{1}{N_{\mathrm{free}}}`$ where $`F_\lambda ^{\mathrm{obs}}`$ is the observed flux density, $`F_\lambda ^{\mathrm{mod}}`$ is the model flux density, $`\sigma ^{\mathrm{obs}}`$ is the rms of the observed flux density. The summation is done over all ULIRG spectra, and $`N_{\mathrm{free}}`$ is the sum of degrees of freedom of modeled spectra. With this indicator the quality of the different models can be compared: the closer to 1., the better. Table 3 presents the mean $`\frac{\chi ^2}{N_{\mathrm{free}}}`$ both for our complete ULIRG sample, and for a high quality subset with observed S/N greater than 3.

#### 4.2.2 Comparison with template sources.

We have obtained model fits for various starburst, AGN, and mixed template sources. A first result can be obtained from a comparison of the two classes of models for the mixed templates Circinus and Centaurus A. Figure ISOCAM-CVF 5-12$`\mu `$m spectroscopy of Ultraluminous Infrared Galaxies 111Based on observations with ISO, an ESA project with instruments funded by ESA Member States (especially the PI countries: France, Germany, the Netherlands and the United Kingdom) with the participation of ISAS and NASA shows the obvious superiority of models where the extinction towards the AGN is different from that towards the starburst. Fully equivalent good fits are achieved for models 3a-3c.

Figures ISOCAM-CVF 5-12$`\mu `$m spectroscopy of Ultraluminous Infrared Galaxies 111Based on observations with ISO, an ESA project with instruments funded by ESA Member States (especially the PI countries: France, Germany, the Netherlands and the United Kingdom) with the participation of ISAS and NASA and ISOCAM-CVF 5-12$`\mu `$m spectroscopy of Ultraluminous Infrared Galaxies 111Based on observations with ISO, an ESA project with instruments funded by ESA Member States (especially the PI countries: France, Germany, the Netherlands and the United Kingdom) with the participation of ISAS and NASA show the observed spectrum and the fitted spectrum computed with model 3a, for several starburst and AGN templates. The Galactic center has been grouped with the AGNs since the very center shown here is also a continum source devoid of UIBs (Lutz et al., 1996). The absence of these features in the GC-SgrA\* spectrum is due to the presence of an H II-region with an intense radiation field that probably destroys the UIB carriers. We overlay the contributions of the starburst component and the contributions of the continuum component. As expected, the starburst component dominates totally the modelled spectrum for the starburst galaxies. Similarly, a quasi-pure continuum contribution is needed to fit the AGN spectra, with extinction being negligible in many cases but significant in others (NGC 5506). The low metallicity galaxy NGC 5253 (not shown) is well fitted by a strong continuum contribution and weak UIBs. The fits to the starburst templates (Fig. ISOCAM-CVF 5-12$`\mu `$m spectroscopy of Ultraluminous Infrared Galaxies 111Based on observations with ISO, an ESA project with instruments funded by ESA Member States (especially the PI countries: France, Germany, the Netherlands and the United Kingdom) with the participation of ISAS and NASA) show at their long wavelength end a varying contribution of a faint steeply rising continuum. Most likely, this is due to variations in the strength of the 15$`\mu `$m VSG (Very Small Grains) continuum, a conclusion also reached by Sturm et al. (2000) on the basis of the SWS spectra of M82 and NGC 253.

A further test is to apply the model to objects which are well known to be composite, containing both starburst and an AGN. In Figure ISOCAM-CVF 5-12$`\mu `$m spectroscopy of Ultraluminous Infrared Galaxies 111Based on observations with ISO, an ESA project with instruments funded by ESA Member States (especially the PI countries: France, Germany, the Netherlands and the United Kingdom) with the participation of ISAS and NASA spectra of four such objects are presented, showing plausible decompositions into the starburst and continuum components. A stringent test can be obtained for Centaurus A and Circinus which have both been observed by ISOCAM/CVF at high spatial resolution (Mirabel et al. 1999; Moorwood et al. 1999 and in preparation). These observations allow us to separate the contribution of the nucleus and the surrounding star-forming regions, and thus a direct test of our model decomposition. The spectrum of the Centaurus A nucleus, extracted from the ISOCAM-CVF observation, is overplotted in Figure ISOCAM-CVF 5-12$`\mu `$m spectroscopy of Ultraluminous Infrared Galaxies 111Based on observations with ISO, an ESA project with instruments funded by ESA Member States (especially the PI countries: France, Germany, the Netherlands and the United Kingdom) with the participation of ISAS and NASA. The modelled continuum contribution is indeed very similar to the observed spectrum of the nucleus, confirming the validity of our approach. We performed the same comparison for the Circinus galaxy (A.F.M. Moorwood, private communication) and find the same good agreement between the observed nuclear continuum and the modelled continuum contribution. At least for these two sources, the separation into two contributions (continuum and starburst) is reliable.

### 4.3 Results of the modelling of starburst and continuum contributions

Table 3 shows that the two classes of models do not reproduce the ULIRG spectra equally well. It is obvious that the models of the second class, where extinctions of continuum and starburst are different, are in better agreement with the observations. This better accuracy is not explained trivially by the increasing number of free parameters. Model 2 which has the same number of free parameters as the ‘good’ models 3a-c is clearly less successful. That means, more flexibility in the shape of the extinction does not improve the fit significantly, while invoking different extinction towards the starburst and AGN does. This lends support to the plausible idea that two different extinctions are applicable to the physically distinct starburst and AGN parts of a ULIRG. Comparing the continuum and starburst extinction derived by our fit procedure, we find the continuum extinction to be slightly smaller on average. This is plausible if the obscuring dust is not distributed with spherical symmetry. We note, however, that this result should be considered tentative given the considerable uncertainty of individual Av values.

Figure ISOCAM-CVF 5-12$`\mu `$m spectroscopy of Ultraluminous Infrared Galaxies 111Based on observations with ISO, an ESA project with instruments funded by ESA Member States (especially the PI countries: France, Germany, the Netherlands and the United Kingdom) with the participation of ISAS and NASA compares the L/C ratio used in previous papers with the ratio of the integrated 5-10$`\mu `$m fluxes contained in the continuum contribution and in the starburst contribution as deduced from the model 3a. In this plot we also add the starburst and AGN comparison samples as described by Rigopoulou et al. (1999). The two indicators are well (anti)correlated. This correlation shows a posteriori the consistency of the two indicators. Table 4 displays the individual value of L/C and the ratio of the integrated 5-10$`\mu `$m fluxes contained in the continuum contribution and in the starburst contribution.

Genzel et al. (1998) and subsequent papers adopted a L/C threshold of 1 to discriminate between starburst and AGN dominated systems. At this threshold, the flux ratio of the continuum and starburst components, each integrated over the 5-10$`\mu `$m range, will be larger than one, simply because the UIBs are relatively narrow features. This is reflected in Figure ISOCAM-CVF 5-12$`\mu `$m spectroscopy of Ultraluminous Infrared Galaxies 111Based on observations with ISO, an ESA project with instruments funded by ESA Member States (especially the PI countries: France, Germany, the Netherlands and the United Kingdom) with the participation of ISAS and NASA. The purpose of the threshold, however, is to identify the component dominating the bolometric luminosity and not just the 5-10$`\mu `$m flux. The mid-infrared contributions have to be extrapolated to the bolometric contributions. Using the spectral energy distributions of M82 and of the central region of NGC 1068 as starburst and AGN templates (Sturm et al. 2000 for M82, including IRAS points and correcting for aperture effects), we find the ratio $`L_{510\mu m}/L_{\mathrm{IR}}`$ about 2.5 times greater in the case of the AGN than for starbursts. Very similar values are obtained for other AGN templates (NGC 4151, 3C273), consistent with the notion of stronger mid-IR continua in Seyferts than in starbursts (de Grijp et al., 1985). Orientation related variations in mid-infrared Seyfert continua (Clavel et al., 2000) caution that there will be some scatter in these bolometric corrections for individual AGN. In the remainder of this paper, we will apply the factor of 2.5 wherever total infrared (i.e. bolometric) quantities are discussed. Taking into account the correction factor 2.5, Figure ISOCAM-CVF 5-12$`\mu `$m spectroscopy of Ultraluminous Infrared Galaxies 111Based on observations with ISO, an ESA project with instruments funded by ESA Member States (especially the PI countries: France, Germany, the Netherlands and the United Kingdom) with the participation of ISAS and NASA confirms the L/C threshold of 1 as a sensible one for differentiating between sources dominated by the starburst component and those dominated by the AGN component.

## 5 Discussion

### 5.1 Mid-infrared spectra and power sources of the highest luminosity ULIRGs

Determination of the power source of the highest luminosity ULIRGs is a prime goal of our project, extending the ISOPHOT-S sample which is dominated by lower luminosity sources (Lutz et al., 1998; Rigopoulou et al., 1999). Figure ISOCAM-CVF 5-12$`\mu `$m spectroscopy of Ultraluminous Infrared Galaxies 111Based on observations with ISO, an ESA project with instruments funded by ESA Member States (especially the PI countries: France, Germany, the Netherlands and the United Kingdom) with the participation of ISAS and NASA displays the L/C ratio for the CAM-CVF ULIRGs together with data from Lutz et al. (1998) and the two hyperluminous sources observed by ISO (Taniguchi et al., 1997; Aussel et al., 1998). The UIB L/C is used to classify sources as starburst-like if L/C $`>`$1 or AGN-like if L/C $`<`$ 1. From our enlarged database, we confirm that most ULIRGs are starburst dominated in the mid-infrared but that there is a trend towards AGN-like ULIRGs at higher luminosities. In our sample which contains 17 sources above $`L_{\mathrm{IR}}=10^{12.4}L_{}`$, the highest luminosity source classified as starburst is at $`10^{12.65}L_{}`$ (IRAS 03835-6432). Such a luminosity limit for the ‘most luminous starburst’ is in agreement with the timescale argument given by Heckman (1994). Adopting the bolometric luminosities of our sources and sizes like those derived by Soifer et al. (2000) for the mid-IR emitting region of ULIRGs, the surface brightness will exceed the limit for starbursts suggested by Meurer et al. (1997). Individual star forming region are known to exceed this limit, however, and the huge gas concentration in the centers ULIRGs provides sufficient fuel. Additionally, Meurer et al. (1997) derived their limit using, for similar sources, H$`\alpha `$-based sizes which are larger than the mid-IR ones found by Soifer et al. (2000), thus staying below the actual surface brightness. Even in our sample enriched in luminous ULIRGs, statistics at the highest luminosities is limited, thus maintaining the possibility for higher luminosity starbursts. Given the modest sample size, it is difficult to assign an accurate value for the luminosity at which half of the sources are starburst and half AGN dominated. From Figure ISOCAM-CVF 5-12$`\mu `$m spectroscopy of Ultraluminous Infrared Galaxies 111Based on observations with ISO, an ESA project with instruments funded by ESA Member States (especially the PI countries: France, Germany, the Netherlands and the United Kingdom) with the participation of ISAS and NASA, this occurs most likely at luminosities $`10^{12.4}L_{}`$ to $`10^{12.5}L_{}`$.

The trend from starburst dominance for most ‘normal’ ULIRGs to AGN dominance at the highest luminosities now consistently emerges from a variety of indicators: Mid-infared continuum (Sanders et al. 1988a, ), mid-infrared spectroscopy (Lutz et al. 1998 and this paper), and optical spectroscopy (Veilleux et al., 1999). As noted by (Lutz et al., 1999), the good agreement between optical and mid-infrared suggests that AGN in ULIRGs do not remain fully embedded for long but manage to break the obscuring screen at least in certain directions.

The highest luminosity ‘starburst-like’ systems appear well above our cutoff L/C of 1, but below the L/C$`3`$ typical for starbursts and lower luminosity ULIRGs. This can be interpreted in different ways. Either they contain a noticeable though not dominant AGN component, or the increased continuum contribution is due to an increase in the H II-region related continuum, as might be reasonable for sources with very high star formation rates in a small region.

### 5.2 The nature of the mid-infrared continuum

The simple L/C method will not be able to break the degeneracy for the origin of ULIRG mid-infrared continua. Also, the parameterization chosen for our fits (Sect. 4.1) describes the spectrum of a ULIRG in terms of just two components : a UIB-dominated starburst component and a continuum component which is either due to an AGN or extra H II region contributions. If the relative weight of ‘PDR-like’ and ‘H II region-like’ components varies for different starburst conditions, as seen e.g. for different regions of the Antennae galaxies (Vigroux et al., 1996), assigning the continuum entirely to an AGN may be misleading. The nature of the continuum has to be investigated in more detail to quantify the fraction of energy radiated by the AGN (AGN continuum) and by star formation activity (UIB-dominated component plus H II region continuum).

It is possible to break the continuum degeneracy by using the slope of the continuum (Laurent et al. 1999). The ratio of the fluxes at 6$`\mu `$m and 15 $`\mu `$m is different for an H II region and a continuum emitted by the dust (torus) surrounding an AGN. Dust in the inner parts of H II regions emits a mid-infrared spectrum that is steep over our wavelength range. This can be seen both in Galactic H II regions (Verstraete et al., 1996; Cesarsky et al., 1996) and in extragalactic regions of particularly concentrated star formation, e.g. the star forming knot A in the Antennae galaxies (Vigroux et al., 1996; Mirabel et al., 1998). In contrast, the continuum emitted by hot AGN dust is much flatter in the mid-infrared, as shown by the examples of Cen A (Mirabel et al., 1999) and NGC1068 (Lutz et al., 2000).

Because of the presence of the UIB emission and silicate absorption features, the steepness of the continuum can be determined reliably only with sufficient wavelength coverage, ideally from 6$`\mu `$m or less to 15 $`\mu `$m (rest wavelength). This condition is not met by the ISOPHOT-S spectra which form most of the ISO database of low resolution ULIRG spectra but cut at less than 11$`\mu `$m rest wavelength for most ULIRGs. It is partially met for the ISOCAM-CVF spectra presented in this paper which typically extend to rest wavelength about 12$`\mu `$m. It is instructive to discuss the case of Arp 220 for which a full CVF spectrum is available (Charmandaris et al., 1997). Arp 220 shows 6.2 and 7.7$`\mu `$m UIBs on no or weak continuum, whereas the continuum in the 11-15$`\mu `$m region is clearly stronger than in a standard M82-type starburst spectrum. This suggests a strong but steep H II continuum as in the star forming knot A of the Antennae, subject to even higher extinction (see Figure ISOCAM-CVF 5-12$`\mu `$m spectroscopy of Ultraluminous Infrared Galaxies 111Based on observations with ISO, an ESA project with instruments funded by ESA Member States (especially the PI countries: France, Germany, the Netherlands and the United Kingdom) with the participation of ISAS and NASA).

The case is less clear for our spectra which do not extend to 15$`\mu `$m. Some spectra suggest a strong but steep continuum at the long wavelength end of a UIB dominated spectrum (e.g. IRAS 03521+0028) while others are fit by a UIB dominated spectrum plus a flat continuum (e.g. IRAS 18030+0705). At the limited wavelength coverage, the uncertainties are considerable, also taking into account that fits are affected by some 7-8$`\mu `$m features not showing the typical PAH shape (Sect. 3).

### 5.3 Model results: Contributions of starburst and AGN activity

Our model can be used to better quantify the contribution of starburst and AGN activity independent from the L/C criterion used earlier. We measured the mid infrared fluxes emitted by starburst and continuum contribution (not dereddened) integrated over the wavelength range 5-10$`\mu `$m (similar to the filter LW2 of ISOCAM). When discussing total infrared (bolometric) contributions instead of 5-10$`\mu `$m ones, the different AGN and starburst SEDs are considered by increasing the starburst contribution by the correction factor 2.5 derived in Sect. 4.3. Since we were not able to reliably break the degeneracy of AGN and H II region continuum for our wavelength coverage and signal-to-noise ratio, we explored two extreme cases for the origin of the continuum, reality most likely being in between.

The first scenario assumes that the continuum is coming entirely from an AGN. The examples of H II region continua discussed above, however, clearly argue for some star formation contribution to the continuum. It is unlikely that an almost pure continuum spectrum is produced in this way. The only star-formation powered galaxies approaching this case are low metallicity dwarfs like NGC 5253, an unlikely match to the ULIRG situation. We hence have adopted the following second extreme scenario: All faint continua, defined as emitting less than 50% of the total 5-10$`\mu `$m flux, are assumed to be due to star formation, while all stronger ones are due to an AGN.

Figure ISOCAM-CVF 5-12$`\mu `$m spectroscopy of Ultraluminous Infrared Galaxies 111Based on observations with ISO, an ESA project with instruments funded by ESA Member States (especially the PI countries: France, Germany, the Netherlands and the United Kingdom) with the participation of ISAS and NASA shows the median ratio of the energy emitted by the starburst component to the total energy emitted ($`F_{510\mu m}(\mathrm{Starburst})`$/$`F_{510\mu m}(\mathrm{Total})`$) of ULIRGs presented in Table 4 and two HYLIRGs (Taniguchi et al., 1997; Aussel et al., 1998), versus the total infrared luminosity. A trend towards AGN dominance at high luminosity is again clearly seen. The cutoff between starburst-dominated galaxies and AGN-dominated galaxies lies at $`L10^{12.412.5}L_{\mathrm{}}`$. The significance of this trend can be tested by dividing the luminosity range into the bins $`L_{\mathrm{FIR}}=10^{11.812.4}L_{}`$ and $`L_{\mathrm{FIR}}=10^{12.413.0}L_{}`$ and using respectively the Kolmogorov-Smirnov’s test and the Student’s t-test to test if the two samples follow the same distribution or have the mean value. The probability for such hypotheses are $`1.0\times 10^4`$(resp. $`1.3\times 10^4`$), demonstrating that the trend is statistically very significant. We also perform the two statistical test while removing the two most luminous galaxies ($`L_{\mathrm{FIR}}>10^{13}L_{}`$). Probabilities are then $`9.0\times 10^4`$ and $`5.0\times 10^4`$.

To quantify the AGN contribution at different luminosities, we identified either all or none of the continuum of each source as AGN according to the two scenarios described previously. We then determined the ratio $`F_{510\mu m}(\mathrm{Starformation})`$/F$`{}_{510\mu m}{}^{}(\mathrm{Total})`$, with $`F_{510\mu m}(\mathrm{starformation})`$ being the sum of UIB-dominated starburst flux and that part of the continuum flux that is ascribed to star formation under the given scenario. For each scenario, we averaged the results for the ULIRGs in each of the two luminosity bins. Table 5 displays the resulting average values along with their uncertainties. The average contribution of star formation to the 5-10$`\mu `$m flux is 65 to 84% for the low luminosity bin and 24 to 33% for the high luminosity bin, the ranges representing the two extreme scenarios. Extrapolating from mid-infrared fluxes to total far-infrared ones, the average contributions of star formation are 82 to 94% for the low luminosity bin and 44 to 55% for the high luminosity bin. These results from quantitative modeling considerably reinforce the conclusions from the earlier work based on the L/C ratio: Low luminosity ULIRGs are predominantly starburst powered, while AGN powering dominates above $`10^{12.4}`$ to $`10^{12.5}L_{\mathrm{}}`$.

## 6 Conclusions

We have presented ISOCAM-CVF low resolution mid-infrared spectroscopy of 16 ultraluminous infrared galaxies obtained within a program targeted at investigating the most luminous ULIRGs. The data are analyzed along with the complete ISO database of low resolution spectra of ULIRGs/HYLIRGs, totalling 76 sources. We use the presence of the mid-infrared ‘UIB’ aromatic bands as a diagnostic of starburst activity, finding starburst-dominated systems up to a luminosity of $`L_{\mathrm{IR}}=10^{12.65}L`$. Other spectra show a strong AGN continuum with weaker features of uncertain origin. We have found one highly obscured AGN.

The mid-infrared spectra can be modeled in terms of a superposition of a UIB-dominated starburst spectrum and a continuum, both potentially strongly obscured. This continuum contains both AGN and additional H II region contributions. The fits prefer a two zone model in which the extinction to starburst and continuum regions may differ. Results from these fits agree well with previous simpler diagnostics based on the line-to-continuum of the UIB features.

Low luminosity ULIRGs are mostly starburst dominated, but the AGN fraction increases with luminosity and most ULIRGs above $`10^{12.4}`$ to $`10^{12.5}L_{\mathrm{}}`$ appear AGN like. We have separated ULIRGs into two luminosity bins of $`L_{\mathrm{FIR}}=10^{11.812.4}L_{}`$ and $`L_{\mathrm{FIR}}=10^{12.413.0}L_{}`$. We find that the average contributions of star formation to the infrared luminosity are 82 to 94% for the low luminosity bin and 44 to 55% for the high luminosity bin.

We are grateful to Michael Rowan-Robinson for allowing us access to the QDOT ULIRGs prior to publication. The work of QDT is supported by Verbundforschung (50 OR 9913 7). DCH acknowledges support from NASA grant NAG5-3359 to US ISO Observers. DR acknowledges support of the EC TMR Network ”Infrared Surveys” (under contract number ERBFMRX-CT96-0068). YT is supported by the Ministry of Education, Science, Sports and Culture of Japan under grants 10044052, and 10304013. SWS and the ISO Spectrometer Data Center at MPE are supported by DLR (DARA) under grants 50 QI 8610 8 and 50 QI 9402 3.

## Appendix A Reduction scheme for CVF observations of faint point sources

Aiming for a good CAM-CVF spectrum including a reliable continuum, a main source of uncertainty is due to the Point Spread Function (PSF). Because the width of the PSF changes with the observed wavelength, it is crucial to determine the observed total flux for each wavelength separately. This is easily done for bright point sources because it is easy to sum the observed flux over a large box (typically 5 by 5 pixels) or to fit the PSF profile to the observed image. Since our sources are point-like at CAM spatial resolution, it is not necessary to consider the complex situation for extended sources. Even for point sources, however, the PSF correction is difficult if they are faint. Then, only a few pixels around the brightest one have a significant S/N (hereafter the brightest pixel will be called “central pixel”). In such a case, the S/N of the corrected spectrum obtained by summing over a box or fitting the PSF profile is usually very low.

A better S/N can be obtained just using the signal of the central pixel and applying a PSF correction factor that varies slowly with wavelength. The main difficulty in determining this factor is to obtain the position of the source within the central pixel. This could be done by fitting a model of the PSF to the observed images at each wavelength. But, because the S/N is quite low for the pixels around the central pixel, direct application of this method gives uncertain positions and a final spectrum where the S/N is low, comparable to the S/N obtained by summing over a box. In order to improve the S/N we adopted the following method :

In CIA, theoretical PSFs are available at every half micron of the ISOCAM-CVF wavelength range, residual variations over such ranges are small. We divided the observed wavelength range into sub-ranges of one micron centered around the wavelength of a theoretical PSF, and summed images corresponding to wavelengths in this range of one micron around the given wavelength (i.e. about 10 images). We then fitted a PSF profile to the resulting image (hereafter the summed image). This method improves the S/N, and the deduced position of the source in the central pixel is known with a better confidence level. We used this position and the PSF to derive the correction factor from flux in the ’central’ pixel to total flux, which was applied for each image in the considered wavelength range.

Although this method gives a better position of the source, it is still ncessary to compute the systematic error introduced by this correction. The formal way is to compute the S/N for each pixel of the summed image, and then to use an inverse Monte-Carlo method to compute the S/N of the spectrum. Another equivalent way which avoids a Monte-Carlo simulation is to consider the images used to obtain the summed image as guesses of a Monte-Carlo simulation . This is equivalent to a simulation with ten guesses because the S/N of the summed image is deduced from these images. We then used the positions of the source obtained for all these images to correct the flux of the summed image. ¿From the scatter of these measurements, it is then possible to deduce the systematic error introduced by the fitting method of the PSF profile.


format: text

text:
# Compton scattering from the proton at NLO in the chiral expansion

## I Introduction

Heavy baryon chiral perturbation theory (HBCPT) is a systematic framework in which to describe the interactions of nucleons with pions and photons. Respecting gauge and Lorentz invariance, it reproduces low-energy theorems, which have been known since the 50’s and 60’s, as the first terms in chiral expansions, and in a number of cases also makes further parameter-free predictions. One of the simplest and cleanest processes to which it can be applied, theoretically at least, is Compton scattering, and indeed this was first studied in the very early days of the theory. First spin-independent scattering was considered. At lowest order (technically second order in HBCPT) the scattering amplitude is simply the energy-independent Thomson term, which is entirely independent of all details of the proton’s structure. This arises simply as a seagull term in the Lagrangian, one of the terms that arise when antinucleons are integrated out of the lowest-order relativistic Lagrangian. Energy and angle dependence in the amplitude enter at third order, where the only diagrams which contribute are one-pion loop graphs. The next lowest terms in an energy expansion of the amplitude (of order $`\omega ^2`$) have coefficient which are called the polarisabilities, $`\alpha `$ and $`\beta `$, and at this order in HBCPT they are simply predicted in terms of well-known nucleon and pion properties . The numerical values, $`\alpha =12.5\times 10^4`$ fm3 and $`\beta =\alpha /10`$, are in remarkably good agreement with the current world averages for the proton, $`\alpha =12.1\pm 0.8\pm 0.5\times 10^4`$ fm3 and $`\alpha +\beta =14.1\pm 0.5\times 10^4`$ fm3 . This was an early triumph for HBCPT.

The calculation of the forward spin-independent scattering amplitude to fourth order followed quickly . However at that order there is no further predictive power for the polarisabilities within the strict framework of the theory, as there are four unknown low-energy constants in the fourth-order Lagrangian which give contributions to the isoscalar and isovector $`\alpha `$ and $`\beta `$. Also, the amplitude is not easily compared with experiment on its own.

The next step was the calculation of the full spin-dependent amplitude to third order. There are four independent amplitudes, each with a leading term (of order $`\omega `$) which is given by a low-energy theorem . The next terms, of $`𝒪(\omega ^3)`$, have coefficients, denoted $`\gamma _i`$, which are called spin polarisabilities. Unfortunately these are not well known experimentally; indeed till recently only the contribution which enters forward scattering, called $`\gamma _0`$, was known at all, and that only from pion photoproduction data rather than direct Compton scattering experiments. More recently $`\gamma _\pi `$ (for backward scattering) has been extracted from unpolarised Compton scattering experiments, with conflicting results . The first direct measurements of polarised photon-proton scattering at MAMI yield a more reliable value $`\gamma _0`$ of $`0.8\pm 0.1\times 10^4`$ fm4  which is within the range of previous estimates. In addition pion photoproduction data has been used as input into fixed-$`t`$ dispersion relations to obtain estimates of the $`\gamma _i`$ individually . Unfortunately the lowest order HBCPT values are in total disagreement with the experimental numbers, the only saving point of agreement being that the isovector quantities, which vanish in LO HBCPT, do seem to be smaller in general than the isoscalar ones.

Most recently, three groups have calculated the spin polarisabilities to fourth order (NLO) . Unlike the spin-independent case, this order contains no LEC’s, and so there are again genuine predictions. Surprisingly the NLO contributions are as large as the LO ones. The agreement with the experimental numbers however is variable.

Of course, classic HBCPT has no explicit $`\mathrm{\Delta }`$. This is assumed to have been integrated out, and its legacy is in the LEC’s of the nucleon-only theory. It is seen first at different orders in different processes: second order for pion scattering, fourth order for spin-independent and fifth order for spin dependent Compton scattering. It is expected to be important; it may improve the agreement of the spin polarisabilities, but risks spoiling the spin-independent ones.

All this emphasis on the polarisabilities obscures the fact the HBCPT does in fact give complete amplitudes for scattering, which one might hope would be believable up to energies of the order of the pion mass, since the chiral expansion is both in powers of $`\omega /\mathrm{\Lambda }`$ and $`M_\pi /\mathrm{\Lambda }`$ (where $`\mathrm{\Lambda }`$ is either $`m_N`$ or the chiral scale $`4\pi f_\pi `$). When the full third-order amplitudes were published, it was claimed that they gave good agreement with cross section data for unpolarised scattering up to the photoproduction threshold. However a more wide-ranging study was done by Babusci et al., who showed that the agreement started to breakdown for backward angles above about 120 MeV. They stopped their comparisons at the threshold.

We have now calculated the full scattering amplitude to fourth order in HBCPT. Since the fourth-order contributions to the spin polarisabilities are so large, it is clearly of interest to see whether the full amplitudes and cross sections are also changed dramatically, and this is the subject of this paper. We have compared with the available data up to 200 MeV , and find that in fact the difference between third- and fourth-order cross sections is small. Above the threshold the good agreement at forward angles continues, while the lack of agreement at backward angles is dramatically worse.

## II The cross section in HBCPT

The usual notation for the scattering amplitude in the Breit frame is, for incoming real photons of energy $`\omega `$ and momentum $`𝐪`$ and outgoing real photons of the same energy and momentum $`𝐪^{}`$,

$`T`$ $`=`$ $`ϵ^\mu \mathrm{\Theta }_{\mu \nu }ϵ^\nu `$ (1)
$`=`$ $`\mathit{ϵ}^{}\mathit{ϵ}A_1(\omega ,\theta )+\mathit{ϵ}^{}\widehat{𝐪}\mathit{ϵ}\widehat{𝐪}^{}A_2(\omega ,\theta )`$ (5)
$`+i𝝈(\mathit{ϵ}^{}\times \mathit{ϵ})A_3(\omega ,\theta )+i𝝈(\widehat{𝐪}^{}\times \widehat{𝐪})\mathit{ϵ}^{}\mathit{ϵ}A_4(\omega ,\theta )`$
$`+\left(i𝝈(\mathit{ϵ}^{}\times \widehat{𝐪})\mathit{ϵ}\widehat{𝐪}^{}i𝝈(\mathit{ϵ}\times \widehat{𝐪}^{})\mathit{ϵ}^{}\widehat{𝐪}\right)A_5(\omega ,\theta )`$
$`+\left(i𝝈(\mathit{ϵ}^{}\times \widehat{𝐪}^{})\mathit{ϵ}\widehat{𝐪}^{}i𝝈(\mathit{ϵ}\times \widehat{𝐪})\mathit{ϵ}^{}\widehat{𝐪}\right)A_6(\omega ,\theta ),`$
where hats indicate unit vectors. By crossing symmetry the functions $`A_i`$ are even in $`\omega `$ for $`i=1,2`$ and odd for $`i=36`$. The amplitudes may be decomposed into pole and non-pole pieces. The pole, which occurs at unphysical values of $`\omega `$ and $`t(=2\omega ^2(\mathrm{cos}\theta 1))`$, arises from Born diagrams in which the intermediate nucleon is on shell, and the contribution containing this pole can be calculated using Dirac nucleons and on-shell photon-nucleon couplings. The leading terms are given by low-energy theorems. With the pole terms truncated at fourth order (that is, at $`𝒪(1/m_N^3)`$), the amplitudes read

$`A_1(\omega ,\theta )`$ $`=`$ $`{\displaystyle \frac{Q^2e^2}{m_N}}+{\displaystyle \frac{e^2}{4m_N^3}}\left((Q+\kappa )^2(1+\mathrm{cos}\theta )Q^2\right)(1\mathrm{cos}\theta )\omega ^2`$ (7)
$`+4\pi (\alpha +\mathrm{cos}\theta \beta )\omega ^2+𝒪(\omega ^4)`$
$`A_2(\omega ,\theta )`$ $`=`$ $`{\displaystyle \frac{e^2}{4m_N^3}}\kappa (2Q+\kappa )\mathrm{cos}\theta \omega ^24\pi \beta \omega ^2+𝒪(\omega ^4)`$ (8)
$`A_3(\omega ,\theta )`$ $`=`$ $`{\displaystyle \frac{e^2\omega }{2m_N^2}}\left(Q(Q+2\kappa )(Q+\kappa )^2\mathrm{cos}\theta \right)+A_3^{\pi ^0}+4\pi \omega ^3(\gamma _1+\gamma _5\mathrm{cos}\theta )+𝒪(\omega ^5)`$ (9)
$`A_4(\omega ,\theta )`$ $`=`$ $`{\displaystyle \frac{e^2\omega }{2m_N^2}}(Q+\kappa )^2+A_4^{\pi ^0}+4\pi \omega ^3\gamma _2+𝒪(\omega ^5)`$ (10)
$`A_5(\omega ,\theta )`$ $`=`$ $`{\displaystyle \frac{e^2\omega }{2m_N^2}}(Q+\kappa )^2+A_5^{\pi ^0}+4\pi \omega ^3\gamma _4+𝒪(\omega ^5)`$ (11)
$`A_6(\omega ,\theta )`$ $`=`$ $`{\displaystyle \frac{e^2\omega }{2m_N^2}}Q(Q+\kappa )+A_6^{\pi ^0}+4\pi \omega ^3\gamma _3+𝒪(\omega ^5)`$ (12)
where the charge of the nucleon is $`Q=(1+\tau _3)/2`$ and its anomalous magnetic moment is $`\kappa =(\kappa _s+\kappa _v\tau _3)/2`$. The polarisability terms are the leading terms in the expansion of the full non-pole amplitudes; the polarisabilities are isospin dependent. The contribution from the $`\pi ^0`$ graph, $`A_i^{\pi ^0}`$, has been separated out, though its contribution is often included in the definition of the polarisabilities. Only four of the spin polarisabilities are independent since three are related by $`\gamma _5+\gamma _2+2\gamma _4=0`$. The spin-dependent amplitudes have no pole contributions beyond the leading ones at this order. (At fifth order there are pieces of order $`\omega ^3/m_N^4`$, though these vanish for both forward and backward scattering amplitudes.) The amplitudes are given for a nucleon spinor normalisation of $`\overline{u}u=1`$.

To calculate these amplitudes in HBCPT, we need the all the diagrams which appear at fourth order. (The pion-loop diagrams are shown in Fig. 1.) There is no clear separation at a diagramatic level between pole and non-pole; heavy-baryon seagulls give pieces of the Dirac pole but, at fourth order, they also contribute to $`\alpha `$ and $`\beta `$. Conversely fourth-order loop diagrams have contributions at $`𝒪(\omega )`$ which correct the bare anomalous magnetic moment in the LET terms. However, taking all diagrams together HBCPT reproduces the pole terms of Eq. 12, with the physical (as opposed to bare) nucleon mass and magnetic moment, as an expansion in powers of $`1/m_N`$. Everything else is non-pole. (While this separation into pole and non-pole is useful in discussion, nothing in this paper depends on them. We simply calculate the full amplitudes to a given order and use them to produce cross sections. This is not affected by the divergence of opinion that is found in the literature over the definition of the polarisabilities .)

For the spin-dependent amplitudes, no LEC’s enter except well-known nucleon and pion properties. We take $`g_A=1.267`$, $`f_\pi =92.32`$, $`\kappa _s=0.120`$, $`\kappa _v=3.706`$ and $`M_\pi =139.6`$ (as the loops are of charged pions in the main). However in the spin-independent amplitudes two sets of LEC’s enter. One set is the pair of LEC’s which contribute to $`\alpha `$ and $`\beta `$; these have been fixed by requiring the full polarisabilities (with third and fourth order contributions) to take the experimental values. The other set are the parameters $`c_i`$ from $`^{(2)}`$, which enter in the tadpole graphs. These are rather poorly known; furthermore they have been determined from pion-nucleon scattering at second, third and fourth order with quite different results (and indeed from different data sets with different results). Clearly the second-order determinations are the most appropriate to use in these calculations as vertex dressing is absent at this order. We have used $`c_1=0.81`$ GeV, $`c_2=2.5`$ GeV and $`c_3=3.8`$ GeV . Sensitivity to these choices is explored later.

All our amplitudes are calculated in the Breit frame. At third order, the result is the same in any frame (so long as the energies remain low); the amplitudes could be in the Breit frame, or the centre-of-mass frame, or even the lab frame (since $`\omega \omega ^{}`$ is of higher chiral order). At fourth order however the result will be different in different frames. Interestingly there is still no dependence on $`\omega \omega ^{}`$, but there will be expressions involving the average nucleon three-momentum $`\overline{𝐩}`$. In the Breit frame $`\overline{𝐩}`$ vanishes, so the amplitudes take on by far their simplest form. In the centre-of-mass frame $`\overline{𝐩}`$ can be written in terms of the photon momenta, and the structures which enter are unchanged, but the amplitudes will have terms with the “wrong” symmetry in $`\omega `$. (These of course can be obtained from a boost of lower order amplitudes, and are not independent.) The differential cross section in any frame is a kinematic prefactor $`\mathrm{\Phi }^2`$ times the invariant $`|T|^2`$, which in terms of the Breit-frame amplitudes is given by

$`|T|^2`$ $`=`$ $`\frac{1}{2}A_1^2(1+\mathrm{cos}^2\theta )+\frac{1}{2}A_3^2(3\mathrm{cos}^2\theta )`$ (16)
$`+\omega ^2\mathrm{sin}^2\theta \left[4A_3A_6+(A_3A_4+2A_3A_5A_1A_2)\mathrm{cos}\theta \right]`$
$`+\omega ^4\mathrm{sin}^2\theta [\frac{1}{2}A_2^2\mathrm{sin}^2\theta +\frac{1}{2}A_4^2(1+\mathrm{cos}^2\theta )+A_5^2(1+2\mathrm{cos}^2\theta )+3A_6^2`$
$`+2A_6(A_4+3A_5)\mathrm{cos}\theta +2A_4A_5\mathrm{cos}^2\theta ]`$
Above the threshold one should replace $`A_i^2`$ by $`|A_i|^2`$ and $`A_iA_j`$ by $`\mathrm{}(A_i^{}A_j)`$.

Experimental data, though taken in the lab frame, are sometimes presented as centre-of-mass frame differential cross sections, so both prefactors are needed:

$$\mathrm{\Phi }_{\mathrm{lab}}=\frac{1}{4\pi }\frac{\omega ^{}}{\omega };\mathrm{\Phi }_{\mathrm{cm}}=\frac{1}{4\pi }\frac{m_N}{\sqrt{s}}.$$
(17)
There is a further complication which stems from the heavy baryon reduction. Physically, there is a threshold at a lab energy of $`\omega _{\mathrm{th}}=M_\pi (1+M_\pi /2m_N)`$. However the HBCPT amplitudes at third order have cusps at $`\omega =M_\pi `$, and at fourth order they are actually singular at this point. The resolution was given by Bernard et al. when discussing forward scattering : if calculated to all orders, the sum of all the singular terms will be precisely what is required to shift the threshold to the physical value. Thus at finite order, to have the threshold at the right point and maintain finite amplitudes, one changes variables from the Breit frame energy $`\omega `$ to $`\zeta \omega /\omega _{\mathrm{th}}`$, then Taylor expands and discards pieces of higher order. At fourth order, therefore, if a naive loop amplitude is $`f^{(3)}(\omega /M_\pi )+f^{(4)}(\omega /M_\pi )`$, the shifted amplitude is

$$f^{(3)}(\zeta )+f^{(4)}(\zeta )+\zeta \left(\frac{\omega _{\mathrm{th}}}{M_\pi }1\right)_{LO}\frac{df^{(3)}(\zeta )}{d\zeta }.$$
(18)
The subscript “LO” indicates the dropping of terms of $`𝒪(1/m_N^2)`$. The shifted amplitude is finite at $`\zeta =1`$.

Since $`|T|^2`$ is an invariant, it should make no difference in which frame it is calculated, and we use the Breit frame as the amplitudes are simplest there. Of course in HBCPT the invariance is only respected up to the order in $`1/m_N`$ to which we are working, and so the results obtained by calculating the amplitudes in different frames would not in fact be identical.

Although we have calculated the amplitudes consistently to fourth order, those amplitudes substituted in Eq. 16 will not give a differential cross section which is of consistent order. Since $`A_1`$ starts at second order with the Thomson term, knowledge of the amplitude to fourth order gives the cross section to sixth order. The fourth order parts of all the other $`A_i`$ are not actually needed, and the seventh- and eighth-order parts of $`|A_1|^2`$ should be subtracted. We do not do this. The various orders of contributions are actually of very similar magnitude, and if we try to remove higher order pieces we may not even be left with a positive cross section. Thus all the cross sections we present are calculated with fourth-order amplitudes, and our kinematic prefactors and change of variable from $`\omega `$ to $`\zeta `$ and from one frame to another are done exactly.

## III Results

In Figs. 2 and 3 we show differential cross sections at various lab and centre-of-mass angles for which data exist below 200 MeV. In each graph we show the Born (pion and nucleon pole) contribution at fourth order, and the HBCPT results with the amplitudes calculated to third and fourth order. The third as well the fourth order results have had the variable change $`\omega /M_\pi \omega /\omega _{\mathrm{th}}`$ as explained above; otherwise the threshold would be in the wrong place. This, and the fact that we take the basic third-order amplitude to be in the Breit frame, rather than the centre-of-mass frame, accounts for the slight differences between our third-order results and those of Babusci et al.. In Fig. 4 we show three third order curves: amplitudes calculated in the Breit frame and shifted, in the Breit frame and unshifted, and in the CM frame and unshifted as in ref. . The differences between these are higher order.

For scattering angles up to 90o the fit at either third or fourth order is very good except at the highest energies, and much better than the Born terms alone. In the cusp region the fourth order is arguably better than the third order at all angles. The good fit to the data above threshold is a feature which has not been remarked before. However for larger angles both curves start to fall below the data as the threshold is approached, and fall dramatically short above. Of course one does not expect agreement to continue indefinitely, as in experiments the $`\mathrm{\Delta }`$ resonance is clearly seen, and this will not be reproduced by HBCPT. However it is not obvious why backward angles should be so much worse than forward angles. As the data in that region is only from one experiment  one might question its accuracy, but in fact it fits well with dispersion analyses which incorporate the more numerous higher energy data.

At fifth order, the new diagrams which enter fall into several groups. First, there are two-loop graphs, which one might well expect to be small. Then, there are one-loop graphs with two insertions from $`^{(2)}`$. Some of these give higher $`1/m_N`$ corrections and these again ought to be small (some indeed have already been incorporated in the change of variables which corrects the threshold) but some, such as two insertions of the anomalous isovector magnetic moment, may be numerically significant. Then there are one-loop graphs with one insertion from $`^{(3)}`$, and seagulls contributing directly to the polarisabilities. It is in these that the effects of the $`\mathrm{\Delta }`$ will show up in the spin-dependent amplitudes, and these contributions are expected to be anomalously large. We have not calculated fifth-order loop graphs. However the seagulls may be incorporated effortlessly, since the free LEC’s will just be fit to give the “experimental” values of the $`\gamma _i`$. In Fig. 5 we show the effect of setting the $`\gamma _i`$ to the values found in three fixed-$`t`$ dispersion-relation studies . Below threshold the agreement is slightly improved, but above the enhancement at backward angles is not sufficient to give agreement.

Finally, we show the effects of varying the input parameters within a purely fourth-order calculation. These are the places where the effect of the $`\mathrm{\Delta }`$ is already present; for instance the tadpole graphs of Fig. 1n-r could arise from integrating out an intermediate $`\mathrm{\Delta }`$. These graphs vanish for forward scattering but contribute significantly to the spin-independent amplitudes at backward angles. First we show the results of varying $`\alpha \beta `$ while keeping the sum fixed. This shows up only at backward angles; to reproduce the Mainz data at 180o we need to make a dramatic shift, to $`\alpha =8`$ and $`\beta =6\times 10^4`$ fm3. (Rather similar numbers were found by Grießhammer and Rupak in a study of Compton scattering from the deuteron in the framework of effective field theory without dynamical pions .) However even these extreme values do not cure the problem above the threshold. Then we show the effect of using a very different parameter set for the $`c_i`$ that enter the tadpole graphs, while keeping the polarisabilities fixed. (These graphs do of course contribute to the polarisabilities, as given ref. . However the fourth order LEC’s are adjusted to give the measured polarisabilities whatever the other contributions might be.) These parameters are quite uncertain and differ substantially in different fits; we show two second-order fits and one fourth-order . Again, the only effect is at backward angles, and none of the parameter sets gives agreement above the threshold.

Thus the conclusion of this fourth-order study is similar to that of the third-order study of Babusci et al. . For angles up to 90o , HBCPT without an explicit $`\mathrm{\Delta }`$ does an excellent job of fitting the experimental differential Compton scattering cross section up to about 170 MeV. However it fails at larger angles for photon energies above about 120 MeV. Very large fifth-order contributions would be needed to approach the experimental data above the threshold at these angles.

This paper was partly supported by the UK EPSRC. I would like to acknowledge the hospitality of the Institute for Nuclear Theory, Seattle, where a substantial part of the work was carried out, and to thank Mike Birse and Daniel Phillips for useful conversations.

## A Fourth-order loop amplitude in the Breit frame

The full amplitudes in the Breit frame for the diagrams of Fig. 1 are as follows. The notation $`𝐭_i`$ is used for the tensor structures which multiply the amplitudes $`A_i`$ of Eq. 5; for example $`𝐭_1=\mathit{ϵ}^{}\mathit{ϵ}`$.

$`T_a`$ $`=`$ $`{\displaystyle \frac{g^2e^2}{4m_Nf_\pi ^2}}(𝐭_1+𝐭_3)[(d1){\displaystyle \frac{J_2[\omega ,m^2]}{\omega }}(w^2+\frac{1}{2}t){\displaystyle \frac{J_0[\omega ,m^2]}{\omega }}]+\text{crossed}`$ (A1)
$`T_b`$ $`=`$ $`{\displaystyle \frac{g^2e^2}{2m_Nf_\pi ^2}}{\displaystyle _0^1}dx[(𝐭_1+𝐭_3)(d+1){\displaystyle \frac{J_2[x\omega ,m^2]}{x\omega }}`$ (A3)
$`(x\omega ^2(𝐭_2+\frac{1}{2}𝐭_5)+(x^2w^2+\frac{1}{2}xt)(𝐭_1+𝐭_3)){\displaystyle \frac{J_0[x\omega ,m^2]}{x\omega }}]+\text{crossed}`$
$`T_c`$ $`=`$ $`{\displaystyle \frac{g^2e^2}{2m_Nf_\pi ^2}}\tau _3(𝐭_1+𝐭_3)\left(\omega J_0[\omega ,m^2]+\mathrm{\Delta }_\pi [m^2]\right)+\text{crossed}`$ (A4)
$`T_d`$ $`=`$ $`{\displaystyle \frac{g^2e^2}{2m_Nf_\pi ^2}}\tau _3(𝐭_1+𝐭_3){\displaystyle _0^1}𝑑x\left(\omega J_0[x\omega ,m^2]+\mathrm{\Delta }_\pi [m^2]\right)+\text{crossed}`$ (A5)
$`T_e`$ $`=`$ $`{\displaystyle \frac{g^2e^2}{2m_Nf_\pi ^2}}(1\tau _3)(𝐭_1+𝐭_3){\displaystyle \frac{1}{\omega }}\left(J_2[\omega ,m^2]J_2[0,m^2]\right)+\text{crossed}`$ (A6)
$`T_f`$ $`=`$ $`{\displaystyle \frac{g^2e^2}{4m_Nf_\pi ^2}}{\displaystyle _0^1}dx[(\mu _v\mu _s\tau _3)(((12x)𝐭_3+𝐭_1)\mathrm{cos}\theta 𝐭_2(12x)𝐭_4)\omega J_0[x\omega ,m^2]`$ (A8)
$`(1\tau _3)(2(d+1)𝐭_1\omega ^1J_2[x\omega ,m^2]+(2x(1x)𝐭_1\frac{1}{2}(12x)𝐭_6)\omega J_0[x\omega ,m^2]+\text{crossed}`$
$`T_g`$ $`=`$ $`{\displaystyle \frac{g^2e^2}{4m_Nf_\pi ^2}}\left((\mu _v+\mu _s\tau _3)((𝐭_1+𝐭_3)\mathrm{cos}\theta 𝐭_2+𝐭_4𝐭_5)+\frac{1}{2}(1+\tau _3)𝐭_6\right)\omega {\displaystyle _0^1}𝑑xJ_0[x\omega ,m^2]`$ (A10)
$`+\text{crossed}`$
$`T_h`$ $`=`$ $`{\displaystyle \frac{g^2e^2}{4m_Nf_\pi ^2}}{\displaystyle _0^1}dy{\displaystyle _0^{1y}}dx[(d+1)(d+3)𝐭_1{\displaystyle \frac{J_2[\stackrel{~}{\omega },m^2xyt]}{\stackrel{~}{\omega }}}`$ (A16)
$`+(𝐭_1(2(d+3)V(x,y)t)+\omega ^2𝐭_2((d+3)((d+5)xy2x2y)+2)+`$
$`((d+3)x1)\omega ^2(𝐭_6𝐭_5)+(d+3)(1xy)\omega ^2𝐭_4){\displaystyle \frac{J_0[\stackrel{~}{\omega },m^2xyt]}{\stackrel{~}{\omega }}}`$
$`+(𝐭_1V(x,y)^2(x𝐭_6x𝐭_5+(1xy)𝐭_4)V(x,y)\omega ^2`$
$`((d+5)xyxy)(2V(x,y)\omega ^2𝐭_2+(1xy)\omega ^4𝐭_7)+xyt\omega ^2𝐭_2)2{\displaystyle \frac{J_0^{}[\stackrel{~}{\omega },m^2xyt]}{\stackrel{~}{\omega }}}`$
$`+(𝐭_2xyV(x,y)^2+𝐭_7xy(1xy)\omega ^2V(x,y))4\omega ^2{\displaystyle \frac{J_0^{\prime \prime }[\stackrel{~}{\omega },m^2xyt]}{\stackrel{~}{\omega }}}]+\text{crossed}`$
$`T_i`$ $`=`$ $`{\displaystyle \frac{g^2e^2}{4m_Nf_\pi ^2}}𝐭_1{\displaystyle _0^1}[(d^21)\mathrm{\Delta }_\pi [m^2x(1x)t]+2(12(d+1)x(1x))t\mathrm{\Delta }_\pi ^{}[m^2x(1x)t]`$ (A18)
$`+4x^2(1x)^2t^2\mathrm{\Delta }_\pi ^{\prime \prime }[m^2x(1x)t]]`$
$`T_j`$ $`=`$ $`{\displaystyle \frac{g^2e^2}{m_Nf_\pi ^2}}{\displaystyle _0^1}dy{\displaystyle _0^{1y}}dx[(d+1)𝐭_1\mathrm{\Delta }_\pi [m^2xyt]`$ (A21)
$`+\left((2V(x,y)t)𝐭_1+2(x+y(d+3)xy)\omega ^2𝐭_2\right)\mathrm{\Delta }_\pi ^{}[m^2xyt]`$
$`+(4V(x,y)2t)xy\omega ^2𝐭_2\mathrm{\Delta }_\pi ^{\prime \prime }[m^2xyt]]`$
$`T_k`$ $`=`$ $`{\displaystyle \frac{g^2e^2}{2m_Nf_\pi ^2}}𝐭_1{\displaystyle _0^1}𝑑x\left[(d1)\mathrm{\Delta }_\pi [m^2x(1x)t]+(12x(1x))t\mathrm{\Delta }_\pi ^{}[m^2x(1x)t]\right]`$ (A22)
$`T_l`$ $`=`$ $`{\displaystyle \frac{g^2e^2}{m_Nf_\pi ^2}}𝐭_1\mathrm{\Delta }_\pi [m^2]`$ (A23)
$`T_m`$ $`=`$ $`{\displaystyle \frac{g^2e^2}{8m_Nf_\pi ^2}}𝐭_1(3\tau _3)(d1)\mathrm{\Delta }_\pi [m^2]`$ (A24)
$`T_n`$ $`=`$ $`{\displaystyle \frac{e^2}{2f_\pi ^2}}𝐭_1\left(8c_3\frac{1}{m_N}(1+\tau _3)\right)\mathrm{\Delta }_\pi [m^2]`$ (A25)
$`T_p`$ $`=`$ $`{\displaystyle \frac{e^2}{2f_\pi ^2}}𝐭_1\left(16c_3\frac{1}{m_N}(1+\tau _3)\right)\mathrm{\Delta }_\pi [m^2]`$ (A26)
$`T_q`$ $`=`$ $`{\displaystyle \frac{8e^2}{f_\pi ^2}}{\displaystyle _0^1}dy{\displaystyle _0^{1y}}dx[\frac{1}{2}(\stackrel{~}{c}_2+(d+2)c_3)𝐭_1\mathrm{\Delta }_\pi [m^2xyt]`$ (A30)
$`+((2c_1M_\pi ^2(\stackrel{~}{c}_2+c_3)\stackrel{~}{\omega }^2+V(x,y)c_3)𝐭_1+`$
$`(xy\stackrel{~}{c}_2+((d+4)xyxy)c_3)\omega ^2𝐭_2)\mathrm{\Delta }_\pi ^{}[m^2xyt]`$
$`+xy\omega ^2𝐭_2(4c_1M_\pi ^2+2\stackrel{~}{\omega }^2(\stackrel{~}{c}_2+c_3)2V(x,y)c_3)\mathrm{\Delta }_\pi ^{\prime \prime }[m^2xyt]]`$
$`T_r`$ $`=`$ $`{\displaystyle \frac{2e^2}{f_\pi ^2}}𝐭_1{\displaystyle _0^1}dx[(dc_3+\stackrel{~}{c}_2)\mathrm{\Delta }_\pi [m^2x(1x)t]`$ (A32)
$`(4c_1M_\pi ^2+2x(1x)tc_3)\mathrm{\Delta }_\pi ^{}[m^2x(1x)t]]`$
$`T_s`$ $`=`$ $`{\displaystyle \frac{3g^2e^2}{8m_Nf_\pi ^2}}(1+\tau _3)(d1)𝐭_1\mathrm{\Delta }_\pi [m^2]`$ (A33)
where the $`c_i`$ are low energy constants from $`_{N\pi }^{(2)}`$, and $`\stackrel{~}{c}_2=c_2g^2/8m_N`$. The integrals $`J_0[\omega ,m^2]`$, $`J_2[\omega ,m^2]`$ and $`\mathrm{\Delta }_\pi [m^2]`$ have their usual meanings, prime denotes differentiation with respect to $`m^2`$, $`\stackrel{~}{\omega }=(1xy)\omega `$ and

$$V(x,y)=\stackrel{~}{\omega }^2+\frac{1}{2}t(1xy+2xy).$$
(A35)
We have also introduced an extra tensor structure, $`𝐭_7=𝝈(\widehat{𝐪}^{}\times \widehat{𝐪})\mathit{ϵ}^{}\widehat{𝐪}\mathit{ϵ}\widehat{𝐪}^{}.`$ This is not independent; $`𝐭_7=\mathrm{sin}^2\theta 𝐭_3+\mathrm{cos}\theta 𝐭_5𝐭_6`$, but as it arises naturally in the calculations—and enters at too high an order in $`\omega `$ to affect the polarisabilities—it is a useful notation.

The notation “+ crossed” means that to every term, another is added with $`\mathit{ϵ}\mathit{ϵ}^{}`$, $`𝐪𝐪^{}`$ and $`\omega \omega `$. Since the $`𝐭_i`$ are all either symmetric or antisymmetric under this transformation, the net effect is to add a term with $`\omega \omega `$ to the coefficients of $`𝐭_1`$ and $`𝐭_2`$, and subtract such a term from the coefficients of $`𝐭_3\mathrm{}𝐭_7`$.

The spin-dependent terms do not all look identical to the expressions given in ref. , but they are equivalent.