text:
# Infrared L Band Observations of the Trapezium Cluster: A Census of Circumstellar Disks and Candidate Protostars

## 1 Introduction

The discovery of small amplitude Doppler periodicities in the spectra of a number of nearby, sunlike, stars has established the existence of extra-solar planets beyond doubt (Marcy & Butler, 1999). However, we do not yet know how common planetary systems are in the Galaxy. Obtaining a meaningful census via direct and indirect detections of planetary systems is a formidable challenge which will likely require the development of new technological capabilities such as precision photometric and astrometric observations from space. On the other hand, an estimate of the frequency of planetary systems is indirectly possible via observations of protoplanetary systems. Such observations are feasible with existing technology, even from the ground.

The low inclinations and eccentricities of the planets around the sun have long suggested that the solar system formed in a disk of dust and gas which surrounded the sun in its earliest stages of evolution. Moreover, theory suggests that circumstellar disk formation is the natural consequence of star formation in rotating collapsing protostellar clouds (Shu, Adams & Lizano, 1987) and that protoplanetary disks should be relatively common in the Galaxy. A protoplanetary disk is considerably easier to detect than a planetary system with a similar mass of solid material (Beckwith & Sargent, 1996). Indeed, circumstellar disks have been directly and indirectly observed around numerous young stellar objects via optical, infrared, centimeter, and millimeter-wave observations (e.g., Beckwith, 1999; Menard & Bertout, 1999; Lada, 1999a; Wilner et al., 2000) and even resolved in optical and millimeter-wavelength images (McCaughrean & O‘Dell, 1996; Dutrey et al., 1996). Thus obtaining a statistically meaningful census of circumstellar disks around young stars may be a very efficient first step toward determining the the frequency of planetary systems in the Galaxy.

One particularly efficient method of producing meaningful measurements of disk frequencies around young stars is to obtain multi-wavelength, infrared imaging surveys of young stellar clusters. Circumstellar disks are bright infrared emitters and produce emission in excess of that emitted by the photospheres of young stars. Such infrared excess is readily measured with near-infrared imaging observations. Young clusters contain hundreds of young target stars spanning nearly the full range of stellar mass within a relatively small angular extent on the sky. Infrared imaging surveys can thus obtain a complete sampling of cluster members in relatively modest observing times.

Infrared color-color diagrams, constructed from multi-wavelength imaging surveys, have been shown to be a useful tool for identifying infrared excesses and candidate circumstellar disks around stars in young clusters (Lada & Adams, 1992). JHK color-color diagrams have been the most commonly employed to investigate infrared excess and circumstellar disk frequency in young clusters. Such studies have produced evidence which suggests that the frequency of infrared (JHK) excess sources decreases with cluster age on timescales on the order of only a few million years, and this potentially sets important constraints on both the lifetimes of circumstellar disks and the duration of the planet building phase within them (Lada & Lada, 1995; Lada, Alves & Lada, 1996; Lada, 1999b). However, JHK observations do not extend to long enough wavelength to enable either a complete or unambiguous census of circumstellar disks in young clusters for two reasons. First, the magnitude of the near-infrared excess produced by a star/disk system depends on the parameters such as the disk inclination, accretion rate and presence and size of inner disk holes (Adams, Lada, & Shu, 1987; Lada & Adams, 1992; Meyer, Calvet & Hillenbrand, 1997), and not all the star/disk systems in a young stellar population are expected to exhibit large enough excesses to be detected in a JHK color-color diagram (Lada & Adams, 1992; Kenyon & Hartmann, 1995; Hillenbrand et al., 1998). Second, a JHK excess can be produced by other phenomena such as extended emission from reflection nebulae and H II regions, which are common in regions of recent star formation (e.g., Haisch, Lada & Lada 2000, hereafter HLL).

The magnitude of the infrared excess produced by a disk rapidly increases with wavelength, and L band (3.5 µm) is likely the optimum wavelength for detecting infrared excesses from protoplanetary disks with ground-based telescopes (HLL). Existing imaging detectors are capable of detecting the photospheres of relatively low mass stars in nearby regions, and circumstellar disks produce sufficiently strong L band excess that almost all stars that have disks can be identified with JHKL (or HKL) observations, independent of inclination, disk accretion rates, etc. (Kenyon & Hartmann, 1995). Moreover, extraneous emission from reflection nebulae, HII regions, etc. is reduced at this wavelength.

This paper is the second in a series which investigate the frequency and evolution of circumstellar disks using L band observations of nearby young clusters. The first paper (HLL) reported results of an L band survey of NGC 2024, a heavily embedded and extremely young cluster in Orion. This paper reports the results of a JHK and L band imaging survey of the famous Trapezium cluster also in Orion.

The Trapezium cluster is a very rich cluster within the Great Orion Nebula. Its members provide a statistically meaningful sample of nearly the entire spread of the stellar mass function. The cluster is somewhat older than NGC 2024 but not as embedded in molecular material, and it is much better studied over a wide range of wavelength (e.g., Herbig & Terndrup, 1986; Churchwell et al., 1987; Felli et al., 1993; McCullough et al., 1995; Stauffer et al., 1994; Ali & DePoy, 1995; Hillenbrand, 1997). There is ample observational evidence for a large population of star/disk systems in the Trapezium cluster. Indeed, the resolved HST images of silhouetted disks within the cluster (McCaughrean & O‘Dell, 1996) are perhaps the most compelling demonstration for the existence of circumstellar disks around young stars yet obtained. However, the size of the disk population in this cluster is uncertain.

Around the O star $`\theta ^1`$C Ori, centimeter-wave radio observations (Churchwell et al., 1987; Felli et al., 1993) and optical imaging both from the ground (Laques & Vidal, 1979) and with HST (O’Dell, Wen & Hu, 1993) revealed what appear to be numerous photo-evaporating disks in the inner region of the cluster. A detailed comparison of these data with deep, high resolution infrared (I & K band) images suggested that the fraction of stars with disks in this part of the cluster is $`>`$ 25-50% and could even be as high as 100% (Stauffer et al., 1994). A detailed analysis of near-infrared IJHK band observations of the larger ONC (Orion Nebula Cluster) by Hillenbrand et al. (1998) suggested a similarly high disk fraction, estimated to be 55-90%. This estimate was supported indirectly by observations of Ca II emission lines which suggested disk accretion occuring in about 70% of the stars in the larger cluster. Infrared L band observations could both confirm and considerably improve upon these estimates of disk frequencies. According to Hillenbrand et al. (1998) only 43 of 1600 visible stars in the ONC cluster have published L measurements and more complete L band observations “…would provide more robust estimates of disk frequencies among both visible and optically obscured young stellar populations.” Earlier, McCaughrean et al. (1996) briefly noted an estimate of 60-80% for the Trapezium disk fraction based on a very rough analysis of L band imaging observations, suggesting that more detailed L band data and analysis could indeed produce a relatively accurate disk census for the cluster. The L band imaging survey reported here confirms these expectations. We derive a significantly improved and very high disk frequency of 80-85 %, for the Trapezium region.

The stars in the Trapezium cluster formed at a relatively high rate (e.g., Lada, Alves & Lada 1996). The cluster is in very close proximity to, but not deeply embedded in, a dense molecular cloud. Its physical association with the molecular cloud raises the interesting issue of the extent to which star formation is still continuing in the cluster region. Recently, Palla & Stahler (1999) have presented evidence indicating that the star formation rate in this cluster has been accelerating with time. However, they postulate a sharp decline in the present day star formation rate due to O star formation. Have most of the stars which will be members of this cluster already formed? Or is rapid star formation continuing to add significant numbers to the cluster population within the remaining molecular gas? Knowledge of the protostellar population in the cloud is essential to address these questions. Protostars are expected to exhibit large L band excesses and very red colors indicative of deeply embedded objects. Our L band observations are both long enough in wavelength and sufficiently sensitive to provide a preliminary reconnaissance of the protostellar population in this cluster/cloud region.

Section 2 discusses the instrumentation employed for this experiment and the data acquisition and reduction procedures used. Section 3 describes the results and analysis of our observations and Section 4 discusses their implications. Section 5 summarizes the primary findings.

## 2 Observations

### 2.1 Instrumentation

Observations of the Trapezium cluster were obtained at the Fred Lawrence Whipple Observatory on Mt. Hopkins, Arizona, using the 1.2m (48in) f/8 Ritchey-Chretien reflector telescope and in conjunction with SAO’s dual channel infrared camera, STELIRCAM (Tollestrup & Willner, 1998). STELIRCAM consists of two 256$`\times `$256 InSb detector arrays. Each is fed from a dichroic mirror that separates wavelengths longer and shorter than 1.9 $`\mu `$m into two independent imaging channels for simultaneous observations on the sky. Three separate magnifications can be selected by moving different cold lens assemblies into the beam. For our JHKL survey we selected a field of view of 2.56′$`\times `$2.56′ with a resolution of 0.6″ per pixel. The magnification is the same for both arrays and both channels simultaneously view approximately the same field on the sky. For all of these observations, the J$`{}_{\text{Barr}}{}^{}(1.25\mu \text{m})`$, H$`{}_{\text{Barr}}{}^{}(1.65\mu \text{m})`$, K$`{}_{\text{Barr}}{}^{}(2.2\mu \text{m})`$, and L$`{}_{\text{Barr}}{}^{}(3.5\mu \text{m})`$ filters were used.

### 2.2 Data Acquisition and Reduction

JHKL band images of the Trapezium cluster were obtained on 14, 15, and 16 December 1997 and JHL images on 04 and 05 November 1998. In both observing runs the cluster was imaged in a square mosaic pattern approximately centered on the O7 star HD 37022 ($`\theta ^1`$C Orionis: R.A.= 05 35 16.47; DEC.=-05 23 22.90; J2000). The individual mosaic images consist of nine separate STELIRCAM fields arranged in a three by three grid. The 1997 fields were spatially overlapped by 34″ ($``$ 25% of the array) in both right ascension and declination providing for accurate positional placement of mosaiced fields as well as redundancy for the photometric measurements of sources located in the overlapped regions. The resulting mosaiced field covered an area of $``$36 arcmin2 (i.e., 6.55′ by 6.55′ ) but because the two arrays were not exactly coincident on the sky the actual area covered simultaneously at all wavelengths was slightly smaller. Typically, the dual channel camera was configured to obtain either simultaneous H and K band or J and L band observations.

Nine off-cluster control fields were also observed in 1997 December to determine accurate sky and flat field calibration frames and measure the colors of field stars. The control fields for these observations were centered at R.A.= 05 26 00; DEC.=-06 00 J2000 and were determined to be free of molecular material by inspection of the Palomar Sky Survey Prints and the 100 micron dust opacity maps of Wood, Myers & Daugherty (1994). Our observing procedure was to intersperse observations of control fields with observations of cluster fields. Typically control field observations were obtained after three consecutive positions in the cluster mosaic were imaged.

Each position in the mosaic was observed in a 9 point dither pattern with small, 12″ offsets between dithers. The integration times and co-additions per dither were 15 seconds with 4 co-additions for J, H and K bands, and 0.1 seconds with 400 co-additions at L band, yielding total integration times of 9 minutes for JHK and 6 minutes at L band. Dark frames were observed at the beginning and end of each night, using the same combinations of integration times and co-additions as the observed images.

Elias et al. (1982) infrared standard stars were observed throughout the night at JHKL wavelengths. These included HD 203856, HD 225023, HD 18881, HD 40335, HD 44612, HD 84800, HD 105601. Magnitudes are reported in the instrumental system of STELIRCAM.

Large temporal background variations at L band (due to window dust, instrument/telescope flexture, and detector anomalies), coupled with the relatively long integration times employed and significant extended L band emission from the HII region made it difficult to determine accurate sky levels across the L band mosaic image. As a result the 1997 data proved to be of relatively poor photometric quality and could not be used to construct a mosaic image of sufficient quality to provide a useful image of extended nebulous emission. Therefore additional L band images were obtained in the telescope nodding mode on 4 and 5 November 1998 to improve the sky subtraction and our ability to faithfully image the extended L band emission from the Orion nebula. In this paper we only use this latter data set for the L band measurements. Short integration J and H observations were simultaneously obtained to extend the dynamic range of observations in those bands and improve the photometry of bright sources.

For each field in the 3 x 3 mosaic grid observed in November 1998 we obtained a single exposure and then nodded the telescope to a sky position 450″ west of the cluster where an identical exposure was obtained. After repeating this sequence for all 9 positions in the grid, the center of the grid was shifted by a small dither ($``$ 5-10″ in a random direction) and the entire mosaic was observed again. This was repeated a total of 18 times. For 9 passes the camera was configured to acquire simultaneous J and L band observations and for the remaining 9 passes H and L observations. The integration times for a single position were 1.0 second with 12 coadds at J and H bands and 0.1 seconds with 100 coadds at L band. The total integration time for each field was 3 minutes in the L band. Overlaps of 20″ between positions in this mosaic provided a total area of approximately 7.0′ by 7.0′ at J,H and L bands, or slightly larger than the 1997 observations. These images were observed only near transit, with a range of airmass of only 1.25 to 1.28. However, since no L band standards were obtained for the November observations, the L band photometry was calibrated internally using 34 stars whose L band magnitudes were independently obtained from observations with the NASA Infrared Telescope Facility (IRTF) on Manua Kea Hawaii in February 2000. The resulting uncertainty in the calibrated L photometry was 0.05 magnitudes. For 14 stars with previous L band photometry in the literature the magnitudes agree with the calibrated November photometry to within the 0.05 magnitude calibration uncertainty. Finally, the entire set of November 1998 photometry was directly compared to the more noisy but independently calibrated 1997 December observations (for which simultaneous L band observations of standards were obtained) and found to agree but within a larger uncertainty due to the poor photometric quality of the 1997 data set.

All data were reduced using standard routines in the Image Reduction and Analysis Facility (IRAF) and the Interactive Data Language (IDL). An average dark frame was constructed from the darks taken at the beginning and end of each night’s observations. This dark frame was subtracted from all target observations to yield dark subtracted images. After appropriate linearization, all dark subtracted target frames were then processed by subtracting the appropriate sky frame and dividing by the flat field. For the 1997 observations, sky frames were constructed from the off-cluster control field images by median filtering the nine dithered observations obtained at each off field position. Flat field frames were constructed by normalizing the median sky frames. For the nodded 1998 observations, sky frames were constructed by median combining the nearest 9 sky positions in time to the target observation. Each sky frame was checked to confirm that all stars had been removed by this process. The sky frames were then normalized to produce flat fields for each target frame. The resulting frames were individually inspected to ensure that the best sky subtraction and flat field had been used. Finally all target frames for a given position in the cluster were registered and combined to produce the final image of each field.

Sources were extracted from each of our images using the automated source extractor DAOFIND in IRAF. The detection threshold was set to be 5 $`\sigma `$ above the average sky noise in an image. Each image was then inspected by eye to recover sources missed by DAOFIND. Photometry was then obtained for all the identified sources using PSF fitting routines in the DAOPHOT package (Stetson, 1987) in IRAF. The primary motivation for using PSF fitting photometry was to mitigate the effect of small scale variations in the background “sky” due to the extended emission of the Orion nebula. This was made possible by use of special sky-fitting routines within the DAOPHOT package as briefly described below.

For the JHK images twenty bright, unsaturated and isolated stars were selected in each of the nine individual images for the construction of a master PSF for that mosaic position. For the L band observations a single PSF was constructed for the entire mosaic image from 30 bright, unsaturated and isolated stars within the mosaic field. This was possible because the PSF did not vary across the mosaic due to the manner in which the observations were obtained. (As described earlier, the entire 3 x 3 mosaic was observed in one pass with short (10 seconds) individual integration times. The final mosaic was constructed from median filtering 18 separate passes.) Sky-fitting routines were employed to improve the subtraction of nebulous emission from the stellar profiles. These routines iteratively obtain estimates for the local sky by removal of target stars and re-measurement of the residual sky (Parker, 1991). Artificial stars added to a very nebulous K band image showed a factor of two improvement in the photometric errors with the sky-fitting option.

All sources were placed in a single positional grid. Absolute positions for our survey stars were established by comparison with the positions for 45 infrared sources determined by K. Luhman (private communication). The relative positional uncertainty in the derived coordinates was found to be 0.2 arc seconds. The entire data base of our survey observations, including complete photometry and positions of all sources will be published in a separate contribution (Muench et al., 2000).

## 3 Results and Analysis

### 3.1 Images

Figure 1a shows our JHK mosaic image of the Trapezium cluster. This is a false color image with J band coded as blue, H band as green and K band as red. The image is displayed in the standard orientation with North at the top and West on the right. The four famous Trapezium stars are in the center of the image. The distinctive reddish nebulosity northwest of the Trapezium stars marks the location of the well known Kleinmann-Low/IRC2 protostellar region. Figure 1b shows the L band mosaic image constructed from our observations. Extended L band emission is evident throughout the region with the well known bright bar prominent in the southeast. The extended L band nebulosity is due to emission from a prominent dust emission feature at 3.3 microns which originates from small grains (PAHs) which are excited by UV radiation from the OB stars in the cluster (e.g., Sellgren, 1981; Duley & Williams, 1981; Geballe et al., 1989). Comparison of the two images shows a clear displacement of the position of the bar in the L band image relative to the JHK image. The L band emission is located further downstream from the exciting Trapezium stars than the J band emission which dominates the emission in the JHK image. While the J band emission traces the ionized gas of the HII region, the L band emission appears to originate within the photodissociation region at the interface between the HII region and the molecular cloud from which the cluster is emerging (Geballe et al., 1989; Roche, Aitken & Smith, 1989).

Another interesting extended feature apparent in the L band map is a long and thin comet-like structure which appears to emanate from one of the infrared sources approximately 56 arc seconds to the east and 14 arc seconds to the north of $`\theta ^1`$C Ori. This source is located within a narrow dark filament at the bottom edge of the well known “dark bay” which is prominently seen in optical images (for example, see the HST WFPC image of Orion Nebula). In Figure 2 we show a close up view of this feature in both the L and J bands. The narrow dark filament or lane is clearly evident in the J band image, although there is no trace of the cometary structure. The cometary feature appears most clearly in the L band image and is extended at a position angle of $`313^o`$ from the infrared source, nearly orthogonal to the dark filament. The cometary tail extends for roughly 21″ or 9400 AU from the star and is very red, hardly visible in the K or H band images. Although morphologically similar to many of the proplyds observed by HST in this nebula, the L band cometary feature is considerably longer than any of the optical proplyds and moreover does not point back to $`\theta ^1`$ C Ori, the primary ionizing source of the proplyd population. In addition the associated stellar component is identified as a non-nebulous star in the proplyd survey of O’Dell and Wong (1996), and it does not appear as a radio source in the Felli et al. (1993 ) survey. On the other hand, this star has a strong L band infrared excess indicative of a circumstellar disk. This cometary feature could be a jet ejected by the star/disk system at its apex, but its very red color would be still difficult to explain. Our data do not allow us to further constrain the nature of this interesting feature. Narrow-band images (e.g., H2) could reveal shock excited emission within this object and test the interesting proposition that it may be a jet.

Within the region surveyed, we detected above our 10 $`\sigma `$ detection limits 548 L sources (m$`{}_{L}{}^{}<12.0`$), 729 K sources (m$`{}_{K}{}^{}<`$ 16.8), 695 H sources (m$`{}_{H}{}^{}<`$ 17.7), and 600 J sources (m$`{}_{J}{}^{}<`$ 18.5). Of these 10 sigma detections, 595 were simultaneously detected at JHK, 490 simultaneously at JHKL and 46 only at HK. Five sources were detected only at L band. Overall a total of 734 distinct sources were detected.

### 3.2 Color-Color Diagrams: The Infrared Excess Fraction

Figures 3 and 4 show the JHK and JHKL color-color diagrams for those (391) sources with both K and L $`<`$ 12 magnitudes. By selecting sources in this way, our survey only includes those sources whose photospheres are detectable at both K and L bands. This insures a meaningful determination of the excess/disk fraction in this cluster. Almost all the stars with K $`<`$ 12 magnitudes have photometric uncertainties less than 0.1 magnitudes in the J, H and K bands. The typical source in our sample is also characterized by photometric uncertainties less than 0.2 magnitudes at L band. Furthermore, we excluded from our sample any source with an uncertainty greater than 0.25 magnitudes at L band.

The solid curve is the locus of colors corresponding to main sequence stars which range in spectral type from early O to M5 (Bessell & Brett, 1988). The solid straight line extending from the main sequence curve is the locus of colors of classical T Tauri stars (Meyer, Calvet & Hillenbrand, 1997). The three dashed parallel lines define the reddening bands for main sequence stars and classical T Tauri stars (CTTS), respectively. These lines are parallel to the reddening vector given by the Cohen et al. (1981) reddening law which we have adopted for this study.

All stars which fall to the right of the main sequence reddening band possess infrared excess and are circumstellar disk candidates. In the JHK diagram 50% $`\pm `$ 7% of the stars are found to display colors indicative of infrared (K band) excess. In the JHKL diagram, 80% $`\pm `$ 7% of the stars are found to have infrared (L band) excess.

The fraction of stars falling in the excess region of the color-color diagram is sensitive to the adopted reddening law and to a lesser extent to the photometric system used to derive the reddening law. We have determined JHK and JHKL excess fractions for two other reddening laws (obtained in different photometric systems) (Koornneef, 1983; Rieke & Lebofsky, 1985) and find JHK excess fractions of 50% and 49% and JHKL excess fractions of 78% and 75%, respectively. These values are within the quoted (statistical) uncertainties of our determinations using the Cohen et al. (1981) extinction law in the STELIRCAM photometric system.

### 3.3 The K-L Distribution: Candidate Protostars

Many of the reddest and most heavily extincted stars will not be detectable at J or H band and will not appear in the color-color diagrams. Therefore it is useful to examine the distribution of K-L colors of the stars in our survey. Figure 5 displays the frequency distribution of K-L color for the 603 stars in our sample detected in both bands. Fifty-seven of these objects were not detected at either J or H. The distribution of the 603 stars detected in both bands has a broad but prominent peak at K-L $``$ 0.9 magnitudes. This distribution is very similar to that in NGC 2024 derived by HLL. There is no peak near zero color similar to that found by Kenyon & Hartmann (1995) for YSOs in Taurus and characteristic of non-excess or diskless (Class III or WTTS) stars. The shift of the peak to relatively large K-L color is primarily due to infrared excess emission and to a lesser extent to reddening. It would require about 25 magnitudes of visual extinction to produce a K-L color of 1.0 magnitude in a naked, diskless star but examination of the color-color diagrams shows that the vast majority of stars in the cluster have extinctions well under 10 visual magnitudes. Thus, for the Trapezium cluster, the presence of infrared excess emission is the primary cause of this shift to larger values of the K-L color.

In the Taurus population of YSOs, objects with K-L $`>`$ 1.5 magnitudes are almost always protostellar in nature (Kenyon & Hartmann, 1995). In order to compile a catalog of potential protostellar objects in the Trapezium region we consider those sources in our survey with similarly large K-L colors. In the Trapezium we find 94 stars with K-L $`>`$ 1.5. Of these, 59 (63%) are also detected in the J and H bands. Twenty-one of these latter objects possess anomalous colors in the JHKL color-color diagram (see discussion in 3.4 below), likely due to source variability, and/or photometry contaminated by unresolved stellar pairs. These sources are not likely to be protostellar objects. However, 28 of the 38 non-anomalous very red sources detected in all four bands are located to the right of the CTTS reddening band above the termination of the (unreddened) CTTS locus. Sources in this part of the color-color diagram possess extreme excess emission which is often characteristic of protostellar (Class I; see Lada (1999a)) objects. Such sources are frequently embedded in bright infrared reflection nebulae which can shift otherwise heavily reddened sources into this region of color-color diagram by depressing the J-H color relative to the K-L color. Sources in this region therefore can be appropriately considered as candidate protostars. The remaining 10 very red sources detected in all bands are located for the most part high within the CTTS reddening band. The colors of these sources are also consistent with those of protostellar objects, and they can also be considered as reasonable candidates for protostellar status. There are 35 sources with K-L $`>`$ 1.5 that were not detected in either J or H band (and cannot be placed on the color- color diagram). These are also candidate protostars. Although some of these sources could be non-protostellar objects that are more embedded than stars of similar K-L color in Taurus, it would take approximately 40 magnitudes of visual extinction to produce a K-L color of 1.5 magnitudes for a naked star and at least 15 magnitudes for a Class II PMS star with a disk. Finally, there are also 5 objects detected only in L band, and limits on their K-L color significantly exceed 1.5 magnitudes. Altogether we identify 78 candidate protostars from the list of sources with K-L colors in excess of 1.5 magnitudes. The positions and JHKL photometry for these 78 candidate protostars are listed in Table 1.

Figure 6 (a and b) compares the spatial distribution of the 78 protostellar candidates with that of the cluster membership as a whole. Inspection of Figure 6a suggests that the protostellar candidates are more spatially confined within the mapped area than are the general cluster members. Figure 6b compares the surface density distribution of the protostellar candidates with that of all other stars in the cluster. The two distributions are clearly different. The surface density distribution of the general cluster membership is dominated by a single peak, coincident with the Trapezium stars. Away from this peak the distribution falls smoothly off with a slight elongation toward the northwest. In contrast, the surface density distribution of the protostellar candidates is characterized by a more spatially extended and elongated maximum which appears double peaked and is clearly offset from the Trapezium stars (and the peak of the cluster surface density distribution). In addition the overall protostellar surface density distribution appears more elongated than that of the general cluster membership with an orientation in the southwest to northeast direction. Figure 7 compares the distribution of the protostellar candidates with that of millimeter-wave continuum emission from dust in the background molecular cloud (Johnstone & Bally, 1999). The protostellar candidate sources have a distribution which is for the most part, coincident with and very similar in extent to that of dust in the well known molecular ridge at the back edge of the cluster. This spatial correspondence with the molecular cloud suggests that this group of sources is physically related to the molecular cloud and more deeply embedded in it than typical cluster members. This is consistent with the identification of these sources as potential protostars.

Because protostellar objects can also be characterized by K-L colors $`<`$ 1.5 magnitudes and can fall into the CTTS region of the color-color diagram, our estimate of the total protostellar population could be a lower limit. On the other hand, it is not clear what fraction of sources located to the right of the plotted CTTS locus and/or with K-L $`>`$ 1.5 are true protostars. We note that a number of proplyds and a few late type stars fall to the right of the CTTS locus (see section 4.0 below) and they are not likely to be protostellar in nature. Whether or not a significant fraction of the candidate protostars in Table 1 are true protostellar objects remains to be determined. Longer wavelength infrared observations with good spatial resolution are necessary to enable a more robust indentification of protostellar objects in this cluster.

### 3.4 Stars with Anomalous K-L Colors

In the JHKL color-color diagram approximately 27 sources (7% of the sample) have colors which place them in “forbidden” regions either to the left of the main sequence reddening band or below and to the right of the CTTS reddening band. (Figure 8 shows the locations of the various zones the JHKL color-color diagram of the Trapezium cluster.) This is almost three times as many sources as the ten sources which fall into similar regions of the JHK color-color diagram. The typical causes for the presence of stars in these two forbidden regions are large photometric errors and contamination of photometry due to unresolved binaries. Of the 27 anomalous sources, 5 are near our L band detection limit and likely characterized by relatively large photometric L band uncertainties. Most of these sources are in the forbidden zone to the left of the main sequence reddening band. Another 4 objects were found to be blended stars or sources which were found to be multiple objects in the Hillenbrand (1997) database but are unresolved in our images. Thus these two effects account for at least a third of the anomalous sources.

For our data set another factor that could result in a star’s placement in one of the forbidden regions is source variability. This is because the L band observations used in our survey were obtained roughly a year after the JHK observations. To evaluate the possibility of variability as a cause of anomalous stellar colors we compared our JHK observations with those in the compilation of Hillenbrand et al (1998). Of the 256 sources with JHK magnitudes in both our catalog and the Hillenbrand et al. compilation, 25 or 10% had fluxes which differed in all three wavelength bands by more than 0.2 magnitudes. Of these 12 (or 5% of the entire sample of 256 sources) differed by more than 0.5 magnitudes in all three bands. Thus source variability could be responsible for placing as many as 5-10% of the sources in our sample in the anomalous regions of the color-color diagram and this factor alone could account for all the sources with anomalous colors. Indeed, 10 of the 27 sources in the two forbidden zones did vary by more than 0.2 magnitudes in the comparison of our JHK data set with that of Hillenbrand et al. Finally, it is very unlikely that variability can account for the vast majority of infrared sources located in the excess region (CTTS and protostellar) of the JHKL color-color diagram since the expected fraction of variables (5-10%) is small compared to the fraction of infrared excess stars and, moreover, already can account for the all sources in the forbidden regions.

## 4 Discussion

Combining L band photometry with JHK observations has enabled us to determine a more robust estimate of the fraction of stars with infrared excess in the Trapezium than has hitherto been possible. Formally we measure an excess fraction of 80% $`\pm `$ 7%. Examination of our control fields suggests that contamination by background sources with mK $`<`$ 12 is insignificant, especially given the large extinction provided by the Orion Molecular Cloud. Thus, the excess fraction we measure corresponds to the excess fraction of the cluster membership. If the infrared excesses originate in circumstellar disks, then the vast majority of stars in the Trapezium cluster are presently surrounded by and consequently were formed with circumstellar disks. Our measurements, therefore, appear to confirm earlier suggestions of a high disk frequency in this young cluster (Stauffer et al., 1994; Hillenbrand et al., 1998).

The assumption of a disk origin for the measured infrared excesses is bolstered by examination of the observed JHKL colors of sources otherwise known to be star/disk systems both in the Taurus dark clouds and the Trapezium cluster itself. For example, HLL analyzed observations obtained and compiled by Kenyon & Hartmann (1995) of known star/disk (class II) sources in the Taurus clouds and found JHKL and JHK excess fractions of 100% and 69%, respectively. We derive very similar statistics from our JHKL colors of the proplyd population in the Trapezium cluster. Proplyds are small but resolved nebulous sources observed in optical HST images around numerous stars in the Trapezium cluster (O’Dell, Wen & Hu, 1993; O’Dell & Wong, 1996). Theoretical considerations indicate that these nebulous objects are most likely photoevaporating envelopes of circumstellar disks (Johnstone, Hollenbach & Bally, 1998; Bally et al., 1998) indicating that the associated stars have such disks. In our imaging survey we found infrared sources associated with 112 of 139 proplyds listed by (O’Dell & Wong, 1996) for the area covered by our images. Of these 112, 96 were simultaneously detected in the JHKL bands and Figure 9 shows their positions in the JHKL color-color diagram. We find 97% of the proplyd stars to have infrared excess. Furthermore we find that only 70% of the proplyds exhibit infrared excess in the JHK color-color diagram. Together these results indicate that JHKL observations are capable of indentifying most, if not all, the circumstellar disk sources in a stellar population. (The close correspondence of JHK and JHKL excess fractions between class II sources in Taurus and proplyds in Orion confirms that source variability has probably not had a significant effect on our ability to measure excess produced by circumstellar disks in the Trapezium cluster.) The overall distributions of Trapezium excess sources in both the JHK and JHKL infrared color-color diagrams (Figures 3 & 4) are very similar to the corresponding distributions of both the proplyd population in the Trapezium cluster (Figure 9) and class II sources in the Taurus cloud (see Figure 8 of HLL). The overall behavior of the Trapezium excess sources in the JHK and JHKL diagrams is also consistent with the predictions of circumstellar disk models (Meyer, Calvet & Hillenbrand, 1997). Taken together these facts appear to provide compelling evidence that the infrared excess sources in the Trapezium cluster represent a population of circumstellar disk systems and that JHKL observations provide an accurate census of that population.

As mentioned earlier, there is a population of radio sources in the cluster whose emission is also thought to arise in the ionized envelopes of photoevaporated disks (Churchwell et al., 1987; Felli et al., 1993). We detected infrared counterparts for 39 of the 50 VLA sources within the boundaries of our imaging survey. Thirty-three were detected in all four infrared bands. Their locations are shown on the color-color diagram in Figure 9. Most of them are associated with proplyds, independently confirming the identification of the radio sources as photoevaporating disks. Not surprisingly, we find that 87% of the VLA sources have infrared excesses. This fraction is somewhat lower than that for the proplyd population and may indicate that some of the radio sources are associated with diskless stars. Indeed, none of the radio sources without infrared excess is associated with a proplyd.

Published spectral types are available for 285 stars within the region we imaged (Hillenbrand et al., 1998). With knowledge of these spectral types an improved estimate of the disk fraction can be made for this cluster. The vast majority (80%) of these spectroscopically observed stars are characterized by spectral types M3 and earlier. Since we have set the boundary of the reddening band at the color of an M5 star, the excess fraction we have derived for the cluster as a whole (80%) from our JHKL observations likely underestimates the true excess fraction. Figure 10 shows the color-color diagram for all the M stars in our imaged region identified by Hillenbrand (1997). Two different boundaries for the main sequence reddening band are drawn at the locations of an M3 and M5 star, respectively. Stars with spectral types later than M3 are represented as solid circles, while stars with spectral types between M0 and M3 are represented by open symbols. The fraction of M0-M3 stars that fall in the infrared excess region bounded by the (appropriate) M3 reddening band is 87% $`\pm `$ 12% and the fraction of M4-M6 stars found in the reddening band bounded by the M5 boundary is 82% $`\pm `$ 17%. Although the same within the respective uncertainties, these fractions are both formally higher than the 80% determined for the whole population measured from the M5 reddening boundary. Moreover, the fraction of O-M3 stars that fall into the excess region as measured from the M3 boundary, is 84% $`\pm `$ 8%.

Although all the above estimates agree within the uncertainties, the results suggest that the actual infrared excess/disk fraction in the Trapezium cluster is close to 85% and somewhat higher than our original formal estimate of 80%. More significantly, the high excess fraction is independent of the spectral type or mass of stars in the cluster. Table 2 lists the excess/disk fractions as a function of spectral type (counting from the M5 boundary). Circumstellar disks appear to form with equally high frequency (within the uncertainties) around all stars with spectral types of F and later, essentially across the entire span of the stellar initial mass function down to the hydrogen burning limit. This is similar to the situation in NGC 2024 where HLL found the disk fraction to also be independent of stellar mass. The earliest type stars in the Trapezium, however, appear to have a lower disk frequency than the later type stars, although the number of stars in the sample with spectral types A or earlier is small.

The overall disk frequency we measure for the Trapezium is essentially the same as that (86%) derived by HLL from similar L band observations of the NGC 2024 cluster which is also located in Orion at the same distance from the sun. Both clusters are very young, although the members of NGC 2024 are more deeply embedded, and suffer more differential extinction than the Trapezium stars. Taken together the measurement of similarly high disk frequencies in the two clusters suggests that disk formation is a natural consequence of star formation for all stellar masses even in the environment of a dense stellar cluster. Moreover, the timescale for significant inner disk evolution must be longer than the age of the Trapezium cluster, about 106 years. In contrast, JHK observations of NGC 2362 and optical spectroscopic studies of the $`\lambda `$ Ori suggest very low disk frequencies ($`<`$ 5%: Alves et al., 2000; Dolan & Mathieu, 1999) for those clusters. This in turn suggests that significant inner disk evolution does occur on timescales shorter than 4-5 $`\times 10^6`$ yrs, the estimated mean ages of these clusters. Because this important timescale constrains the duration of planet building within circumstellar disks, a more robust determination of it is clearly desirable. We anticipate improved estimates once L band observations of more clusters in our sample are analyzed.

There are 78 potential protostars in the Trapezium out of some 603 stars detected in both K and L bands (Table 1). The fraction of protostellar (candidate) sources in the Trapezium is about 13%, which is similar to the protostellar candidate fraction of about 14% derived from L band observations of NGC 2024 by HLL. For a cluster age of 0.8-1.0 $`\times `$ 106 years, the lifetime of the protostellar stage would be approximately 105 years in the Trapezium if all the candidate sources turned out actually to be protostars. This is typical (within factors of 2) of protostellar lifetimes found in Taurus, Ophiuchus, and NGC 2024. Stated another way, if a significant fraction (i.e., $``$ 50%) of the protostellar candidates in the Trapezium are true protostars, with lifetimes comparable to those in other regions, then our observations are consistent with a relatively constant star formation rate in the cluster up to the current epoch. Thus, star formation could be continuing in the molecular cloud behind the Orion Nebula a rate similar to that which gave rise to the foreground cluster! The molecular ridge behind the Trapezium is part of a more extended molecular cloud which connects with the OMC 1, OMC 2, and OMC 3 clouds to the north. The extended ridge is also a site of intense protostellar activity containing at least a dozen Class 0 objects and outflow sources (Chini et al., 1997; Yu, Bally & Devine, 1997; Johnstone & Bally, 1999). This suggests a similarly high rate of star formation in the extended ridge, given that the lifetimes of Class 0 objects are about an order of magnitude shorter than those of Class I protostars. Recently Palla & Stahler (1999) have re-analyzed the ages for the stars in the Orion Nebula Cluster and have suggested that the rate of star formation in this region has been accelerating with time with the most vigorous episode taking place in the last two million years. They have also suggested that the formation of O stars in the cluster should squelch further star formation with the result that the star formation rate should be drastically decreased at present or in the near-future. Palla and Stahler estimate that the timescale for this steep decline should be less than 106 years. If the above estimate of the protostellar population is close to the true size, then it appears that the predicted rapid decline has not yet occured (primarily because the O stars have not completely cleared the region of dense star forming gas). However, our identifications of protostars in the Trapezium region are not secure, and longer wavelength observations are necessary to obtain a more accurate estimate of the protostellar population in the cluster and a more robust determination of the current star formation rate in the cloud behind it.

## 5 Summary

We have obtained a sensitive L band (3.4 µm) survey of a $``$ 36 arcmin2 region centered on the Trapezium cluster in Orion. We have used these observations in conjunction with new JHK imaging observations of the same region to produce a relatively unambiguous census of the circumstellar disk population in the cluster and to compile a list of candidate protostars within the region. The primary results and conclusions derived from our study can be summarized as follows:

1. We detected 391 stars in the surveyed region with K and L $`<`$ 12.0 magnitudes. From analysis of the JHKL color-color diagram we find that 80% $`\pm `$ 7% of the stars exhibit detectable excess emission and are likely surrounded by circumstellar disks. A subsample of 285 of these stars for which there are published spectral types, has a slightly higher excess/disk fraction of 85%. More significantly, within our statistical uncertainties, this high excess/disk fraction is independent of spectral type for stars with spectral types between F-M. Thus the probability of disk formation around a star at the time of its birth is both extremely high and independent of stellar mass across essentially the entire stellar mass spectrum down to the hydrogen burning limit. The disk fraction (42%) for the higher mass O, B and A stars is lower however, and may suggest either a lower probability for disk formation or more rapid disk dispersal times for these stars. Overall the results presented in this paper confirm and improve upon previous studies which also suggested a high disk fraction in the cluster.

2. Of known proplyds in our surveyed region, 97% have infrared L band excess. The fraction of radio continuum emitting stellar systems with L band excess is found to be 87%. These results provide additional confirmation for our identification of infrared excess sources with circumstellar disks, and they illustrate the effectiveness of L band observations for obtaining an accurate census of star/disk systems in the cluster.

3. The disk fraction in the Trapezium cluster is very similar to that (86% $`\pm `$ 8%) recently obtained from a sensitive L band survey of NGC 2024, the other rich and very young embedded cluster in Orion (Haisch, Lada & Lada 2000). Thus our results confirm the findings of the NGC 2024 study, which suggested that, independent of mass, stellar birth in the environment of a rich and dense cluster does not diminish the probability of disk formation around a new star. Indeed, in both these clusters, the vast majority of stars, spanning essentially the entire range of the stellar mass function, were born with circumstellar disks and consequently the potential to form planetary systems. Moreover, the high disk frequencies found in these two clusters indicate that disk lifetimes are on the order of the age of the clusters, $`10^6`$ years, or longer.

4. Seventy-eight of the sources in our survey possess K-L colors typical of deeply embedded (A$`{}_{V}{}^{}2040`$ magnitudes) young stellar objects. Only about half of these stars are detectable at the J and H bands, but most of these are characterized by extreme infrared excess falling outside the CTTS region (see Figure 3 & 4) of the JHKL color-color diagram. The spatial distribution of these very red stars is very similar to that of the dense molecular ridge behind the cluster. These properties are characteristic of protostellar objects, and we identify these objects as candidate protostars (Class I objects). If the size of the protostellar population is even a modest fraction of our candidate protostar sample, then star formation may be continuing within the molecular ridge at a rate similar to that which produced the stars in the foreground Trapezium cluster. Sensitive, longer wavelength observations with good spatial resolution are necessary to determine the fraction of L band candidates that are truly protostellar in nature and evaluate the present level of star forming activity in the Trapezium region.

We thank Doug Johnstone for providing us with the 850 $`\mu `$m map of the OMC 1 region. EAL acknowledges support from a Research Corporation Innovation Award and a Presidential Early Career Award for Scientists and Engineers (NSF AST 9733367). AM is a Smithsonian Predoctoral Fellow. K.E.H. gratefully acknowledges support from a NASA Florida Space Grant Fellowship and an ISO grant through JPL #96104. We also acknowledge support from an ADP (WIRE) grant NAG 5-6751.

Figure Captions


format: text

text:
# Entropy of thermal quasiparticles in nuclei

## I Introduction

Nuclear structure at low excitation energy depends on the residual long-range two-body interaction. One of the most exciting consequences of this interaction is the forming of $`J=0`$ nucleon pairs, the so-called Cooper pairs, where nucleons are moving in time reversed orbitals. It has been a great challenge to observe the breaking of these Cooper pairs and the general quenching of pair correlations as function of temperature $`T`$ and angular frequency $`\omega `$. Rotational high spin states are essentially described by quasiparticles coupled to the rotating underlying core. Backbending phenomena and spin alignments have been described according to their quasiparticle properties.

The concept of thermal quasiparticles is based on the idea that thermal properties at low excitation energy are governed by a few quasiparticles coupled to the cold core of Cooper pairs. These quasiparticles are thermally scattered on available single particle states. Thus, they are not well-defined by one single (Nilsson) orbital with given spin and parity, but exhibit an average of the spectroscopical properties of the orbitals in the vicinity of the Fermi surface. At higher temperatures ($`T>0.5`$ MeV), the pair correlations are quenched and the core of Cooper pairs is no longer well defined .

The entropy is a fundamental quantity in the description of many body systems. For hot nuclei the entropy has the same importance as the spin alignment has for rotational nuclei. It can be evaluated from mechanical variables like energy, particle number and volume. Thus, it can be deduced in the microcanonical ensemble, which is the appropriate ensemble for an isolated system like the nucleus. The microcanonical entropy is defined for all isolated systems, from few-body systems up to infinite systems. It has been shown that phase transitions in small systems can be studied using the microcanonical entropy without invoking the thermodynamical limit .

Although less appropriate, also the canonical ensemble theory is widely applied to nuclei. One advantage of the canonical ensemble is that the thermodynamical observables are smooth functions of temperature due to the Laplace transformation involved in the calculation of the canonical partition function. Therefore, we will model our findings in this work within the canonical ensemble theory.

In a recent paper , we extracted the microcanonical single quasiparticle entropy in hot 161,162Dy and 171,172Yb nuclei to be $`1.70(15)`$ in units of the Boltzmann constant $`k_B`$. Furthermore, the single quasiparticle entropy was found to be approximately constant in the 1–6 MeV excitation region.

In this paper we will investigate the global systematics of nuclear entropy as function of nuclear mass, excitation energy (or temperature), number of excited quasiparticles and ground state spin. In Sect. II we sketch a simple model for thermal quasiparticles in nuclei. The method of determining level density by counting levels and using neutron resonance spacings is described in Sect. III. In Sect. IV we display the systematics of entropies and in Sect. V the properties of thermal quasiparticles are discussed. In Sect. VI we present an application involving thermal quasiparticles and, finally, concluding remarks are given in Sect. VII.

## II Simple model for thermal quasiparticles

In this section we will describe a simple model demonstrating the properties of thermal quasiparticles, which in turn will be compared to experimental data. The main task is to describe the partition function $`Z`$ in the canonical ensemble for a given temperature $`T`$ and number of thermal quasiparticles $`n`$. The thermodynamical quantities of interest can be deduced from the Helmholtz free energy

$$F(T)=T\mathrm{ln}Z(T),$$
(1)
where we express the temperature $`T`$ in units of MeV and $`k_B`$ is set to unity. From this bridge equation connecting statistical mechanics and thermodynamics, we may calculate the entropy $`S`$, the average excitation energy $`E`$, the heat capacity $`C_V`$, and the chemical potential $`\mu `$ at the temperature $`T`$ by

$`S(T)`$ $`=`$ $`\left({\displaystyle \frac{F}{T}}\right)_V`$ (2)
$`E(T)`$ $`=`$ $`F+TS`$ (3)
$`C_V(T)`$ $`=`$ $`\left({\displaystyle \frac{E}{T}}\right)_V`$ (4)
$`\mu (T)`$ $`=`$ $`{\displaystyle \frac{F}{n}},`$ (5)
where $`n`$ is the number of thermal quasiparticles.

In the evaluation of the partition function $`Z`$, we assume spin $`1/2`$ fermions scattered into a doubly degenerated single particle level scheme with equal energy spacing $`ϵ`$. The level scheme has infinitely many levels, and it is one scheme for protons and one for neutrons, as indicated in Fig. 1. These conditions simulate the Nilsson single particle scheme for deformed nuclei.

We will first construct the partition functions $`z_n^{}`$ for $`n`$ identical spin-up fermions placed in one and the same single particle scheme. We assume that $`n`$ quasiparticles are available, with the core of Cooper pairs as a reservoir. The total partition functions can then be expressed by these basic partition functions. The scattering of one spin-up fermion into the level scheme gives the partition function

$$z_1^{}=\underset{i=0}{\overset{\mathrm{}}{}}\mathrm{exp}(iϵ/T)=\frac{1}{1\mathrm{exp}(ϵ/T)}.$$
(6)
When two identical spin-up fermions are placed into the same level scheme, the summing is restricted due to the Pauli principle by

$`z_2^{}`$ $`=`$ $`{\displaystyle \underset{i=0}{\overset{\mathrm{}}{}}}\left[\mathrm{exp}(iϵ/T){\displaystyle \underset{j=i+1}{\overset{\mathrm{}}{}}}\mathrm{exp}(jϵ/T)\right]`$ (7)
$`=`$ $`{\displaystyle \frac{1}{1\mathrm{exp}(2ϵ/T)}}{\displaystyle \frac{\mathrm{exp}(ϵ/T)}{1\mathrm{exp}(ϵ/T)}}`$ (8)
$`=`$ $`z_1^{}{\displaystyle \frac{\mathrm{exp}(ϵ/T)}{1\mathrm{exp}(2ϵ/T)}}.`$ (9)
Generally, for $`n`$ identical spin-up fermions we find

$$z_n^{}=z_{n1}^{}\frac{\mathrm{exp}((n1)ϵ/T)}{1\mathrm{exp}(nϵ/T)}.$$
(10)
We now allow spin up and down for the fermions, and evaluate the corresponding partition function $`z_n`$. For one fermion with spin up or down we simply obtain a degeneration of two, since $`z_1^{}=z_1^{}`$:

$$z_1=z_1^{}+z_1^{}=2z_1^{}.$$
(11)
If two fermions occupy the same level scheme, the partition function becomes more complicated. With the spins of the two fermions antiparallel ($`m=m_1+m_2=0`$), no Pauli blocking is present and we obtain

$$z_2(m=0)=z_1^{}z_1^{}=\left(z_1^{}\right)^2.$$
(12)
The contribution with parallel spin invoke the Pauli principle giving

$$z_2(m=\pm 1)=z_2^{}+z_2^{}=2z_2^{}.$$
(13)
It is now straightforward to evaluate $`z_n`$ for $`n`$ fermions, allowing spin up and down. For convenience we list the seven lowest partition functions here:

$`z_1`$ $`=`$ $`2z_1^{}`$ (14)
$`z_2`$ $`=`$ $`2z_2^{}+\left(z_1^{}\right)^2`$ (15)
$`z_3`$ $`=`$ $`2z_3^{}+2z_2^{}z_1^{}`$ (16)
$`z_4`$ $`=`$ $`2z_4^{}+2z_3^{}z_1^{}+\left(z_2^{}\right)^2`$ (17)
$`z_5`$ $`=`$ $`2z_5^{}+2z_4^{}z_1^{}+2z_3^{}z_2^{}`$ (18)
$`z_6`$ $`=`$ $`2z_6^{}+2z_5^{}z_1^{}+2z_4^{}z_2^{}+\left(z_3^{}\right)^2`$ (19)
$`z_7`$ $`=`$ $`2z_7^{}+2z_6^{}z_1^{}+2z_5^{}z_2^{}+2z_4^{}z_3^{}.`$ (20)
The partition function for nuclei with protons and neutrons may now be constructed from the partition functions $`z_n`$. In Fig. 1 the relevant configurations for the lowest quasiparticle excitations in the proton-odd system are illustrated. For this one-quasiparticle system, there is only one proton partition function to consider. For the three-quasiparticle case the $`\pi \nu ^2`$ and the $`\pi ^3`$ partition functions have to be included. For the two-quasiparticle case of the even-even system we assume that either a $`\pi `$ or $`\nu `$ pair can be broken. In the four-quasiparticle case there are $`\pi ^2\nu ^2`$, $`\pi ^4`$ and $`\nu ^4`$ configurations to take into account. From similar arguments for higher number of quasiparticles, we construct the seven lowest partition functions for proton-odd and even-even systems by

$`Z_1`$ $`=`$ $`z_1`$ (21)
$`Z_2`$ $`=`$ $`2z_2`$ (22)
$`Z_3`$ $`=`$ $`z_3+z_1z_2`$ (23)
$`Z_4`$ $`=`$ $`2z_4+(z_2)^2`$ (24)
$`Z_5`$ $`=`$ $`z_5+z_3z_2+z_1z_4`$ (25)
$`Z_6`$ $`=`$ $`2z_6+2z_4z_2`$ (26)
$`Z_7`$ $`=`$ $`z_7+z_5z_2+z_3z_4+z_1z_6.`$ (27)
The partition functions for the neutron-odd system are identical to the functions for the proton-odd system.

The partition functions for the odd-odd system can be constructed in a similar manner. The two-quasiparticle partition function is simply $`(z_1)^2`$, giving an entropy of twice the single quasiparticle entropy of the odd-mass system. Higher numbers of quasiparticles yield the following partition functions for the odd-odd case:

$`\stackrel{~}{Z}_2`$ $`=`$ $`(z_1)^2`$ (28)
$`\stackrel{~}{Z}_4`$ $`=`$ $`2z_3z_1`$ (29)
$`\stackrel{~}{Z}_6`$ $`=`$ $`2z_5z_1+(z_3)^2`$ (30)
$`\stackrel{~}{Z}_8`$ $`=`$ $`2z_7z_1+2z_5z_3.`$ (31)
The single quasiparticle entropy is easily computed from our model with Eqs. (1,2,27) as

$$S_1(T)=\frac{ϵ/T}{\mathrm{exp}(ϵ/T)1}+\mathrm{ln}\frac{2}{1\mathrm{exp}(ϵ/T)}.$$
(32)
At higher excitation energies more quasiparticles have to be considered and the mathematical expressions become more complicated. In Fig. 2 the calculated $`n`$-quasiparticle entropy $`S_n`$ is shown as function of $`ϵ/T`$. The entropy is seen to scale rather well with the number of quasiparticles, and there is a weak Pauli-blocking effect. The reason for this appearent extensivity is the onset of new degrees of freedom with increasing number of quasiparticles. In our formalism, this is manifested in the additional terms in the partition function appearing from the inclusion of spin up/down quasiparticles and the breaking of proton/neutron pairs.

This approximate extensivity of the entropy can be utilized to estimate the many quasiparticle partition function as a simple product of $`n`$ quasiparticle partition functions $`Z_n=(Z_1)^n`$. The quasiparticle entropy is then roughly given by

$$S_n=nS_1,$$
(33)
where we argue in Sect. IV that the $`S_1`$ value is experimentally available for most nuclei.

To create quasiparticles costs energy. In order to break one Cooper pair of the core the energy $`2\mathrm{\Delta }`$ is necessary. Taking this energy into account we can finally write the partition functions of even-even (ee), odd (odd) and odd-odd (oo) nuclei as weighted sums of the $`n`$-quasiparticle partition functions:

$`Z^{\mathrm{ee}}`$ $`=`$ $`1+e^{2\mathrm{\Delta }/T}Z_2+e^{4\mathrm{\Delta }/T}Z_4+\mathrm{}`$ (34)
$`Z^{\mathrm{odd}}`$ $`=`$ $`Z_1+e^{2\mathrm{\Delta }/T}Z_3+e^{4\mathrm{\Delta }/T}Z_5+\mathrm{}`$ (35)
$`Z^{\mathrm{oo}}`$ $`=`$ $`\stackrel{~}{Z}_2+e^{2\mathrm{\Delta }/T}\stackrel{~}{Z}_4+e^{4\mathrm{\Delta }/T}\stackrel{~}{Z}_6+\mathrm{}`$ (36)
provided an ideal core of Cooper pairs, not changing with the number of quasiparticles. If the system is extensive, we may write the partition functions in closed form. By using the property of extensivity, $`Z_n=(Z_1)^n`$, and performing the sum up to $`N`$ broken Cooper pairsFor $`N\mathrm{}`$ broken Cooper pairs our simple model breaks down since the partition functions become infinite at the temperature $`ϵ/T=\mathrm{ln}[12\mathrm{exp}(\mathrm{\Delta }/T)]`$., we obtain

$`Z^{\mathrm{ee}}`$ $`=`$ $`{\displaystyle \frac{1(e^{\mathrm{\Delta }/T}Z_1)^{2N+2}}{1(e^{\mathrm{\Delta }/T}Z_1)^2}}`$ (37)
$`Z^{\mathrm{odd}}`$ $`=`$ $`Z_1Z^{\mathrm{ee}}`$ (38)
$`Z^{\mathrm{oo}}`$ $`=`$ $`\left(Z_1\right)^2Z^{\mathrm{ee}}.`$ (39)
The latter partition functions are good approximations for the expressions of Eq. (36) for up to $``$2 broken pairs, however, for more than $``$4 quasiparticles, the Pauli blocking cannot be neglected.

Figure 3 shows the entropy as function of temperature calculated from our model, allowing up to 5 broken nucleon pairs. Also calculations assuming extensivity of the entropy, in the sense $`Z_n=(Z_1)^n`$, are shown. There are serious discrepances between the non-extensive and extensive approaches for $`T>0.3`$ MeV.

There are limitations of the present model. The excited states are assumed built by a few quasiparticles, while the other particles are forming Cooper pairs in the underlying core. This picture will break down for (i) light nuclei, (ii) around closed shells and (iii) at high excitation energy. The first two points rest on the fact that the single particle energy spacing should be small or comparable to the pairing interaction in order to build a superfluid phase. The last point is connected to the thermal quenching of pair correlations, where many valence particles participate in the excitation. A fingerprint of the breakdown of the thermal quasiparticle picture is when the entropy in the odd-mass and even-even systems becomes essentially equal. These limitations have been addressed in Ref. .

## III Determination of level density anchor points

The nuclear level density $`\rho `$ is scarcely known as function of excitation energy. Recently, there has been published compilations of level density parameters. In our work we show that it is possible to determine two level density points quite reliably for numerous nuclei. To demonstrate the technique, Fig. 4 shows the two extracted anchor points (filled data points) for 172Yb, together with the level density deduced from known discrete levels (solid lines). The first point $`(E_1,\rho _1)`$ is based on counting known discrete levels at low excitation energy, giving an anchor point in the excitation energy region of $`0<E<2`$ MeV. Information on known levels is taken from the database of Ref. , where we have analyzed levels of $`1600`$ nuclei. The second point $`(E_2,\rho _2)`$ is estimated from the average neutron resonance spacing at the neutron binding energy $`B_n`$, where we use the compilation of Iljinov et al..

Figure 4 also includes the data points (open cirles) obtained with the method of Ref.. The dashed straight line describes the data rather well. This line will be referred to as the constant-temperature level-density formula, defined by

$$\rho (E)=Ce^{E/T}$$
(40)
with $`T^1=(\mathrm{ln}\rho _2\mathrm{ln}\rho _1)/(E_2E_1)`$ and $`C=\rho _1\mathrm{exp}(E_1/T)`$. This line connecting the two anchor points is used to determine the nuclear temperature $`T`$.

In the extraction procedure of $`(E_1,\rho _1)`$, an initial slope is predicted for the straight line; in the case of 172Yb, the inverse of $`T=0.57`$ MeV. Then the two highest values of the count-based $`\mathrm{ln}\rho `$ relative to this line is foundWe pick the highest values since it is more probable that levels have escaped detection than the opposite.. With maximums at the energy bin numbers $`a`$ and $`b`$, we find the lower anchor point by

$$\rho _1=\mathrm{exp}\left[\frac{1}{ba+1}\underset{i=a}{\overset{b}{}}\mathrm{ln}\rho (i)\right]$$
(41)
with energy $`E_1`$ in bin number $`(a+b)/2`$. In the evaluation of the uncertainty of $`\rho _1`$, the number of levels $`N`$ ($`\sigma _N=\sqrt{N}`$) and the uncertainty of $`T`$ ($`\sigma _T=0.1T`$) are taken into account. Nuclei far from $`\beta `$-stability may have been poorly investigated due to low reaction cross sections, and the extracted level density and uncertainties of these nuclei should be adopted with caution.

The second point $`(E_2,\rho _2)`$ is estimated from the average neutron resonance spacing at $`B_n`$. Assuming neutron $`\mathrm{}=0`$ capture on a nucleus with target spin and parity $`I^\pi `$, levels with the same parity and two spins $`I\pm 1/2`$ are populated. In order to transfer the spacing $`\overline{D}`$ between these levels to the level density for all spins and parities, we use the spin distribution of Gilbert and Cameron

$$\rho _2=\frac{2\sigma ^2}{\overline{D}}\left[(I+1)e^{\frac{(I+1)^2}{2\sigma ^2}}+Ie^{\frac{I^2}{2\sigma ^2}}\right]^1.$$
(42)
The spin cut-off parameter $`\sigma `$ is defined through $`\sigma ^2=0.0888A^{2/3}\sqrt{aU}`$, where the back-shifted energy is $`U=EE_{bs}`$. The level density parameter $`a`$ and the back-shift parameter $`E_{bs}`$ are defined by $`a=0.21A^{0.87}`$ MeV-1 and $`E_{bs}=C_{bs}+\mathrm{\Delta }`$, respectively, where the back-shift correction is given by $`C_{bs}=6.6A^{0.32}`$ according to Ref. . The pairing gap parameter is estimated by $`\mathrm{\Delta }=12A^{1/2}`$ MeV.

In Figs. 5-7 we show the anchor points extracted for nine isotopes of germanium, gadolinium and uranium, respectively. The above procedure has been tested carefully and found reliable, as demonstrated in these figures. In Table 1 the two anchor points are listed together with the temperature $`T`$ for the nuclei where both anchor points are known. We also see from Table 1 and Figs. 5-7 that the extracted value of $`T`$ is approximately equal for neighboring nuclei, thus, giving confidence to the methods applied for establishing the two anchor points.

## IV Extraction of entropy

The measured level density $`\rho (E)`$ is proportional to the number of states accessible to the nuclear system at excitation energy $`E`$. Thus, the entropy in the microcanonical ensemble is given by

$$S(E)=\mathrm{ln}(\rho (E)/\rho _0)=\mathrm{ln}\rho (E)+S_0.$$
(43)
Applied to Eq. (40), parameter $`T`$ is interpreted as the thermodynamical temperature in the microcanonical ensemble, since by definition:

$$T(E)=\frac{1}{(S/E)_V}.$$
(44)
The constant shift parameter $`S_0=\mathrm{ln}\rho _0`$ can be adjusted to fulfill the third law of thermodynamics; $`S0`$ when $`T0`$. The condition $`T=0`$ probably holds for the rotational ground band of deformed even-even nuclei. However, also other interpretations for $`T=0`$ are possible, and we will not fix the $`S_0`$ parameter to a specific value in this work, but instead discuss the relative entropy $`SS_0=\mathrm{ln}\rho `$.

The level density in between the two anchor points is not very well known. Recently, the level densities for 161,162Dy, 166Er and 171,172Yb have been measured. A common feature for these nuclei, is that the level densities develop more like the constant temperature level density formula of Eq. (40), than the back-shifted Fermi-gas level density , as pointed out in the discussion of Fig. 4 and Ref. . This constant temperature behavior is interpreted as the breaking of Cooper pairs with increasing excitation energy.

In the following, we will discuss the microcanonical entropy for given excitation energy $`E`$. In order to establish a common play ground, we extrapolate the first and second level density anchor points to excitation energies 1 and 7 MeV, using the constant temperature formula of Eq. (40). Figures 8 and 9 display the entropy $`SS_0`$ based on the two anchor points evaluated at 1 and 7 MeV, respectively. The data are plotted as function of the mass number, where nuclei along the $`\beta `$-stability line are chosen.

## V Thermal quasiparticle properties

### A Regions of quasiparticle extensivity

In the lower panel of figure 8 the difference in entropy between odd-mass and even-even systems $`\mathrm{\Delta }S=S_{\mathrm{odd}}S_{\mathrm{ee}}`$ is shown. The data show fluctuations which are partly due to shell closures, incomplete knowledge of discrete levels at low excitation energy, and structural changes between various isotopes and isotones. For the mid-shell regions the odd-mass nuclei have about $`\mathrm{\Delta }S1.7`$ higher entropy compared to the even-even system. This feature holds for odd-even as well as for even-odd nuclei, however, neutrons carries slightly more entropy than protons. In Fig. 9 the data are sensitive to the level density extracted from neutron resonance level spacing. The fluctuations are less pronounced at 7 MeV because the level density is based on many quasiparticle states, smearing out structural effects. Also here, as indicated on the figure by lines, an entropy difference of $`\mathrm{\Delta }S1.7`$ is found. Here, we should remember that these findings are based on a constant temperature level density extrapolation. If a Fermi gas level density is applied, e.g. the model of Egidy et al. shown in Fig. 4, the extracted microcanonical entropy difference $`\mathrm{\Delta }S`$ will be somewhat modified; it will increase at 1 MeV and decrease at 7 MeV. However, as demonstrated in Ref. , the experimental value of $`\mathrm{\Delta }S`$ is rather constant for excitation energies above 1.5 MeV.

It is amazing that $`\mathrm{\Delta }S`$ is roughly the same both at 1 and 7 MeV of excitation energy. This means that even with 4 – 6 quasiparticles excited at 7 MeV , the last thermal quasiparticle still carries $`\mathrm{\Delta }S1.7`$. Thus, the experimental entropy scales with the number of quasiparticles. This property is also accounted for in our model§§§We should keep in mind that the model is formulated in the canonical ensemble, where a certain temperature $`T`$ corresponds to a wide excitation energy range.. The regions of quasiparticle extensivity can be identified in Figs. 8 and 9 as entropy-gaps between the entropies of odd-mass and even-even nuclei. The extensivity regions are found for mid-shell nuclei with mass number $`A>40`$.

In odd-odd nuclei, several neutron energy spacing data are known . As seen in Fig. 9 (see lines), a proton and a neutron give an entropy excess of around $`2\mathrm{\Delta }S3.4`$, as expected from the discussion in Sect. II. A similar study from known discrete levels is difficult, due to experimental limitations in building level schemes with several hundreds of levels per MeV. However, in the case of the well studied odd-odd 152,154Eu isotopes the lower anchor point seems to have been determined with acceptable confidence. Figure 10 shows that an additional neutron outside the odd-proton core gives an entropy difference as high as $`2.5`$. The extracted temperatures are all lying within a narrow window of $`T=0.570.60`$ MeV.

The feature that $`\mathrm{\Delta }S`$ is approximately independent of excitation energy $`E`$ is shown in Table 1. Here, constant temperature $`T`$ for neighboring mid-shell nuclei means that the slope of $`\mathrm{ln}\rho `$ is equal and, thus, $`\mathrm{\Delta }S`$ is approximately independent of $`E`$. Fig. 11 confirms the systematics, where neighboring nuclei exhibit the same $`(S/E)`$ value. This indicates that the level density, and hence the entropy, has the same functional form as function of excitation energy. However, the feature disappears at shell closures.

### B Constancy of quasiparticle entropy

Figure 12 shows the entropy difference for a neutron outside its even-even core for various rare earth elements. Although various mid-shell isotopes reveal quasiparticle entropy that deviates from $`\mathrm{\Delta }S1.7`$ it is still surprising how constant this quantity appears to be. However, as shown in Fig. 12, the quasiparticle around the $`N=82`$ shell closure carries low entropy due to the limited number of quantuum states (multiplicity) available for the system.

Our model gives a qualitative interpretation of the constancy of the quasiparticle entropy. With the assumptions given, the model shows that the entropy $`S_n`$ is a function of $`ϵ/T`$, only. In order to deduce experimentally the value of $`ϵ/T`$, we recognize that the single particle level spacing $`ϵ`$ is connected to the level density parameter $`a`$ by

$$a=\frac{\pi ^2}{6}(g_p+g_n)\frac{\pi ^2}{3}g=\frac{\pi ^2}{3ϵ},$$
(45)
where the single-particle level-density parameters for protons and neutrons ($`g_p`$ and $`g_n`$) are assumed to be approximately equal. Thus, we obtain using Eq. (44)

$`{\displaystyle \frac{ϵ}{T}}={\displaystyle \frac{\pi ^2}{3a}}\left({\displaystyle \frac{S}{E}}\right).`$ (46)
In Fig. 13 the quantity $`(S/E)/a`$ is plotted as function of mass number $`A`$, using the $`a`$ values of Ref.  and the data of Fig. 11. For increasing $`A`$, the extracted values are found in the range of $`0.100.08`$, which corresponds to $`ϵ/T=0.330.26`$ according to Eq. (46). This variation gives a minor effect on the predicted single particle entropy $`S_1`$; Fig. 2 shows that the variation in $`ϵ/T`$ corresponds to an increase of $`S_1`$ from $`2.7`$ to $`3.0`$.

The experimentally extracted single quasiparticle entropy $`\mathrm{\Delta }S`$ shows no clear systematic trend as function of $`A`$ within the uncertainties. Essentially, the local variations due to structural changes between neighboring nuclei are more pronounced than the weak variation in $`S_1`$ given by the model.

The quasiparticle entropy $`\mathrm{\Delta }S`$ could in principle depend on the odd-mass ground state spin, since an orbital with high spin has more coupling combinations to the underlying rotating core. For excitation energies around the neutron binding energy $`B_n`$, the ground state orbital participates no more than other orbitals around the Fermi level, so no resemblence of the ground state spin is expected. However, around the lower anchor point ($`E_1,\rho _1`$), the ground state orbital could play a role. In Fig. 14 we see no significant effect for rare earth nuclei, and the same result applies to other mass regions. This indicates that many orbitals participate to build up the level density in the lower excitation region of odd-mass nuclei.

To close this subsection, we point out that Fig. 14 indicates an average quasiparticle entropy $`\mathrm{\Delta }S`$ of $`1.5`$ and $`1.9`$ for protons and neutrons, respectively. The valence protons and neutrons in this figure occupy orbitals in the $`Z=5082`$ and $`N=82126`$ shells, containing 32 and 44 orbitals, respectively. Thus, the single particle level energy spacing is expected to scale roughly as $`ϵ_p/ϵ_n=44/321.4`$. According to Fig. 2, a 40% variation in $`ϵ/T`$ around the value $`ϵ/T=0.3`$ gives a change in $`S_1`$ of about 0.4. Therefore, this effect could be the reason why the protons carry slightly less entropy compared to the neutrons.

## VI Application of thermal quasiparticles

There are several applications for the concepts presented here. In analogy with the description of rotational nuclei, the Helmholtz free energy $`F`$ plays the same role as the energy in the rotating frame $`E^{}`$. To visualize the analogy, the two cases read

$`F`$ $`=`$ $`ETS`$ (47)
$`E^{}`$ $`=`$ $`E\omega I_x,`$ (48)
where $`\omega `$ is the rotational frequency and $`I_x`$ is the spin alignment along the rotational x-axis. Generally, $`S`$ and $`I_x`$ are positive quantities, so that $`F`$ and $`E^{}`$ decreases with increasing $`T`$ and $`\omega `$, respectively. A decoupled Cooper pair with high $`I_x`$ may, at a certain critical frequency $`\omega _c`$, become more favoured than its coupled counterpart, giving rise to the rotational backbending phenomenon.

In a similar manner, we will here demonstrate that the free energy of a system with two additional excited quasiparticles becomes favoured at a critical temperature $`T_c`$. This temperature depends on the pairing gap $`\mathrm{\Delta }`$ and the quasiparticle entropy $`\mathrm{\Delta }S`$. If both these quantities would be independent of the number of quasiparticles excited, a break-up of all Cooper pairs should happen at one and the same critical temperature $`T_c`$, giving rise to a strong quenching of the pair correlations.

Let us assume that the partition function for an even-even nucleus is described by $`Z(T)=Z_{\mathrm{core}}(T)`$. This partition function represents the lowest free energy $`F(T)`$ at low temperatures $`T`$. We now introduce the breaking of one pair relative to this coreThe core described by $`Z_{\mathrm{core}}`$ is not necessary a cold core, but may already have pairs broken.. Assuming the two quasiparticles to be independent of the core, we may describe the excited system by a factorization of the partition function into

$$Z^{}Z_{\mathrm{core}}Z_2e^{2\mathrm{\Delta }/T},$$
(49)
where $`Z_2`$ is the two quasiparticle partition function. With this assumption, the difference in Helmholtz free energy becomes

$$\mathrm{\Delta }F=T(\mathrm{ln}Z^{}\mathrm{ln}Z_{\mathrm{core}})=T\mathrm{ln}Z_2+2\mathrm{\Delta },$$
(50)
and factorizing $`Z_2`$ assuming extensivity, gives

$$\mathrm{\Delta }F=T\mathrm{ln}Z_1^2+2\mathrm{\Delta }=2(F_1+\mathrm{\Delta }),$$
(51)
where $`F_1`$ is the free energy of one quasiparticle. This corresponds to the chemical potential defined (see Eq. (5)) as the energy required to bring one quasiparticle out of the reservoir:

$$\mu (T)=\mathrm{\Delta }F/2=F_1+\mathrm{\Delta }.$$
(52)
We have applied this formalism to the experimental data of 161,162Dy, see Ref. . The partition functions for each nucleus are calculated by

$$Z(T)=\underset{E}{}\mathrm{\Delta }E\rho (E)e^{E/T}$$
(53)
where $`\mathrm{\Delta }E`$ is the experimental energy bin. The assumptions made for evaluating $`Z`$ is given elsewhere. The corresponding Helmholtz free energy $`F`$ for 162Dy is displayed in the upper panel of Fig. 15. The upper energy curve, based on data for 161Dy, is given by $`F(^{161}`$Dy) + $`\mathrm{\Delta }_n`$, where the pairing gap parameter is determined from neutron separation energies $`𝒮_n`$ by

$$\mathrm{\Delta }_n=\frac{1}{4}|𝒮_n(N+1,Z)2𝒮_n(N,Z)+𝒮_n(N1,Z)|,$$
(54)
giving $`\mathrm{\Delta }_n=`$ 0.917 MeV for 162Dy. Thus, the energy curve of 161Dy can be interpreted as the free energy for an even-even system with one extra quasiparticle. We therefore identify $`F_1+\mathrm{\Delta }_n`$ in Eq. (51) with the difference $`F(^{161}\mathrm{Dy})+\mathrm{\Delta }_nF(^{162}\mathrm{Dy})`$.

In the lower panel of Fig. 15 the chemical potential $`\mu =F_1+\mathrm{\Delta }_n`$ is shown as function of temperature. The higher entropy in 161Dy compared to 162Dy decreases the free energy relative to the even-even core with increasing temperature. When $`\mu (T)0`$, thermal quasiparticles can be formed without consuming free energy. Thus, from Fig. 15 we estimate the critical temperatureAt the critical temperature $`T_c`$ we have $`F_1=\mathrm{\Delta }=T_c\mathrm{ln}Z_1`$, and applying $`Z_1=2z_1^{}`$, we get $`ϵ/T_c=\mathrm{ln}[12\mathrm{exp}(\mathrm{\Delta }/T_c)]`$, the same expression as in footnote 1. for the quenching of pair correlations to be $`T_c=0.42`$ MeV. The uncertainty of this number is mainly due to the extensivity assumed in the evaluation of $`\mathrm{\Delta }F`$. Furthermore, the adoption of a chemical potential $`\mu =\mathrm{\Delta }`$ at $`T=0`$, rests on the assumption that the Fermi level coincides with single particle energies in both nuclei. Also the extrapolation of the experimental level density to higher excitation energies, see Ref. , gives systematical errors. From these considerations we estimate a 15 % error in the $`T_c`$ value.

The deduced value of $`T_c=0.42(6)`$ is considered to be common to the two dysprosium isotopes analyzed. Recently, values of $`T_c=`$ 0.52(4) and 0.49(4) MeV were found for 161Dy and 162Dy, respectively, using the canonical heat capacity as thermometer for the phase transition.

Finally, we point out a curious connection between the experimental quasiparticle entropy $`\mathrm{\Delta }S`$ and the critical temperature $`T_c`$ deduced from a Fermi gas model with pairing. If we assume a constant $`\mathrm{\Delta }S`$ and a constant differences $`\mathrm{\Delta }`$ between the $`E`$ functions for the two curves of Fig. 15, we find

$$T_c=\frac{1}{\mathrm{\Delta }S}\mathrm{\Delta },$$
(55)
giving $`T_c=0.59\mathrm{\Delta }`$ for $`\mathrm{\Delta }S=1.7`$. This estimate coincides almost exactly with the theoretical expression $`T_c=0.57\mathrm{\Delta }`$ given in Refs. . The theoretical factor $`0.57`$ is independent on the mass number, consistent with the property of the experimental single quasiparticle entropy.

## VII Conclusions

We have analyzed the level densities of nuclei based on counting discrete levels and from neutron scattering data. It was possible to obtain level density anchor points for about 280 nuclei. The level densities have been transformed into microcanonical entropy, revealing the phase space of the system at fixed excitation energies.

The single quasiparticle entropy has been extracted from the odd-even entropy differences. The quasiparticles not coupled in Cooper pairs, exhibit an entropy of $`\mathrm{\Delta }S1.7`$ which is relatively independent of the presence of other quasiparticles. This apparent extensivity is valid for mid-shell nuclei and for all ground state spins. A simple model for thermal quasiparticles in the canonical ensemble also supports the property of extensivity.

The single quasiparticle entropy $`\mathrm{\Delta }S`$ for mid-shell nuclei shows no clear mass dependency. This is explained in our model by the fact that $`\mathrm{\Delta }S`$ only depends on the ratio between the single particle spacing $`ϵ`$ and the temperature $`T`$, which is found to be almost independent of the mass number for mid-shell nuclei.

For thermal systems, the Helmholtz free energy is analogous to the energy in the internal frame of rotating nuclei. This opens for the study of various quantities and critical phenomena as function of temperature. We therefore believe that the concept of single quasiparticle entropy may provide a fruitful basis for the understanding of hot nuclei.

The simple model presented here is only adequate for mid-shell nuclei. Around closed shells one expects exciting effects from the increasing single particle energy spacings. This will also influence the entropy difference between odd-mass and even-even nuclei. Therefore, a statistical description of the transition to closed shells would be of great interest.

We wish to acknowledge the support from the Norwegian Research Council (NFR).


text:
# Dark Matter And Detector Cross Sections

## 1 Calculational Details

In order to get reasonably accurate results, it is necessary to include a number of theoretical corrections in the analysis. We list here the main ones used in the calculations below: (i) In relating the theory at $`M_G`$ to phenomena at the electroweak scale, the two loop gauge and one loop Yukawa renormalization group equations (RGE) are used, iterating to get a consistent SUSY spectrum. (ii) QCD RGE corrections are further included below the SUSY breaking scale for contributions involving light quarks. (iii) A careful analysis of the light Higgs mass $`m_h`$ is necessary (including two loop and pole mass corrections) as the current LEP limits impact sensitively on the relic density analysis for tan$`\beta 5`$. (iv) L-R mixing terms are included in the sfermion (mass)2 matrices since they produce important effects for large tan$`\beta `$ in the third generation. (v) One loop corrections are included to $`m_b`$ and $`m_\tau `$ which are again important for large tan$`\beta `$. (vi) The experimental bounds on the $`bs\gamma `$ decay put significant constraints on the SUSY parameter space and theoretical calculations here include the leading order (LO) and approximate NLO corrections. We have not in the following imposed $`b\tau `$ (or $`tb\tau `$) Yukawa unification or proton decay constraints as these depend sensitively on unknown post-GUT physics. For example, such constraints do not naturally occur in the string models where $`SU(5)`$ (or $`SO(10)`$) gauge symmetry is broken by Wilson lines at $`M_G`$ (even though grand unification of the gauge coupling constants at $`M_G`$ for such string models is still required).

All the above corrections are under theoretical control except for the $`bs\gamma `$ analysis where a full NLO calculations has not been done. (We expect that while the full analysis might modify the regions of parameter space excluded by the $`bs\gamma `$ experimental constraint, the minimum and maximum values of $`\sigma _{\stackrel{~}{\chi }_1^0p}`$ would probably not be significantly changed.) The analysis of $`\sigma _{\stackrel{~}{\chi }_1^0p}`$, taking into account the above theoretical corrections has now been carried out by several groups obtaing results in general agreement . These results are presented below.

Accelerator bounds significantly limit the SUSY parameter space. In the following we assume the LEP bounds  $`m_h>104(100)`$ GeV for tan$`\beta `$=3(5) and $`m_{\chi _1^\pm }>104(100)`$ GeV. (For tan$`\beta >5`$, the $`m_h`$ bounds do not produce a significant constraint.) For $`bs\gamma `$ we assume an allowed range of $`2\sigma `$ from the CLEO data . The Tevatron gives a bound of $`m_{\stackrel{~}{g}}270`$ GeV( for $`m_{\stackrel{~}{q}}m_{\stackrel{~}{g}}`$).

Theory allows one to calculate the $`\stackrel{~}{\chi }_1^0`$-quark cross section and we follow the analysis of  to convert this to $`\stackrel{~}{\chi }_1^0p`$ scattering. For this one needs the $`\pi N`$ $`\sigma `$ term,$`\sigma _{\pi N}`$ and $`\sigma _0=\sigma _{\pi N}(m_u+m_d)p|\overline{s}s|p`$ and the quark mass ratio $`r=m_s/(1/2)(m_u+m_d)`$. We use here $`\sigma _{\pi N}=65`$ MeV, from recent analyses  based on new $`\pi N`$ scattering data, $`\sigma _0=30`$ MeV and r$`=24.4\pm 1.5`$.

## 2 mSUGRA model

We consider first the mSUGRA model where the most complete analysis has been done. mSUGRA depends on four parameters and one sign: $`m_0`$ (universal scalar mass at $`M_G`$), $`m_{1/2}`$ (universal gaugino mass at $`M_G`$), $`A_0`$ (universal cubic soft breaking mass) , tan$`\beta =H_2/H_1`$ (where $`H_{2,1}`$ gives rise to (up, down) quark masses) and $`\mu /|\mu |`$( where $`\mu `$ is the Higgs mixing parameter in the superpotential, $`W_\mu =\mu H_1H_2`$). One conventionally restricts the range of these parameters by “naturalness” conditions and in the following we assume $`m_01`$ TeV, $`m_{1/2}600`$ GeV (corresponding to $`m_{\stackrel{~}{g}}1.5`$ TeV, $`m_{\stackrel{~}{\chi }_1^0}240`$ GeV), $`|A_0/m_0|5`$, and 2$`tan\beta `$ 50. Large tan$`\beta `$ is of interest since SO(10) models imply tan$`\beta 40`$ and also $`\sigma _{\stackrel{~}{\chi }_1^0p}`$ increases with tan$`\beta `$. $`\sigma _{\stackrel{~}{\chi }_1^0p}`$ decreases with $`m_{1/2}`$ for large $`m_{1/2}`$, and thus if one were to increase the bound on $`m_{1/2}`$ to 1 TeV ($`m_{\stackrel{~}{g}}2.5`$ TeV), the cross section would drop by a factor of 2-3.

The maximum $`\sigma _{\stackrel{~}{\chi }_1^0p}`$ arise then for large tan$`\beta `$ and small $`m_{1/2}`$. This can be seen in Fig.1 where ($`\sigma _{\stackrel{~}{\chi }_1^0p}`$)max is plotted vs. $`m_{\stackrel{~}{\chi }_1^0}`$ for tan$`\beta `$=20, 30, 40 and 50. Current detectors obeying Eq (1) are then sampling the parameter space for large tan $`\beta `$, small $`m_{\stackrel{~}{\chi }_1^0}`$ (and also small $`\mathrm{\Omega }_{\stackrel{~}{\chi }_1^0}h^2`$) i.e

$$tan\beta \stackrel{>}{}25;m_{\stackrel{~}{\chi }_1^0}\stackrel{<}{}90\mathrm{GeV};\mathrm{\Omega }_{\stackrel{~}{\chi }_1^0}h^2\stackrel{<}{}0.1.$$
(3)
To discuss the minimum cross section, it is convenient to consider first $`m_{\stackrel{~}{\chi }_1^0}\stackrel{<}{}150`$ GeV ($`m_{1/2}350`$) where no coannihilation occurs. The minimum cross section occurs for small tan$`\beta `$. One finds

$$\sigma _{\stackrel{~}{\chi }_1^0p}\stackrel{>}{}4\times 10^9\mathrm{pb};m_{\stackrel{~}{\chi }_1^0}\stackrel{<}{}140\mathrm{G}\mathrm{e}\mathrm{V}$$
(4)
which would be accessible to detectors that are currently being planned (e.g. GENIUS).

For larger $`m_{\stackrel{~}{\chi }_1^0}`$, i.e. $`m_{1/2}\stackrel{>}{}150`$ the phenomena of coannihilation can occur in the relic density analysis since the light stau, $`\stackrel{~}{\tau }_1`$, (and also $`\stackrel{~}{e}_R`$, $`\stackrel{~}{\mu }_R`$) can become degenerate with the $`\stackrel{~}{\chi }_1^0`$. The relic density constraint can then be satisfied in narrow corridor of $`m_0`$ of width $`\mathrm{\Delta }m_0\stackrel{<}{}25`$ GeV, the value of $`m_0`$ increasing as $`m_{1/2}`$ increases . Since $`m_0`$ and $`m_{1/2}`$ increase as one progresses up the corridor, $`\sigma _{\stackrel{~}{\chi }_1^0p}`$ will generally decrease.

We consider first the case $`\mu >0`$. One finds in general that $`\sigma _{\stackrel{~}{\chi }_1^0p}`$ also decreases as $`A_0`$ increases. Fig.2 shows $`\sigma _{\stackrel{~}{\chi }_1^0p}`$ in the domain of large $`A_0`$ and for two values of tan$`\beta `$. One sees that the smaller tan$`\beta `$ still gives the lower cross section, though the difference is mostly neutralized at larger $`m_{1/2}`$. (For large tan$`\beta `$, $`m_0`$ also becomes large to satisfy the relic density constraint i.e $`m_0700`$ GeV for tan$`\beta `$=40, $`m_{1/2}=600`$ GeV.) We have in general for this regime

$`\sigma _{\stackrel{~}{\chi }_1^0p}`$ $`\stackrel{>}{}`$ $`1\times 10^9\mathrm{pb};\mathrm{for}m_{1/2}600\mathrm{G}\mathrm{e}\mathrm{V},`$
$`\mu `$ $`>0`$ $`,A_04m_{1/2}.`$ (5)
This is still within the sensitivity range of proposed detectors.

When $`\mu `$ is negative an “accidental” cancellation can occur in part of the parameter space in the coannihilation region which can greatly reduce $`\sigma _{\stackrel{~}{\chi }_1^0p}`$ . This can be seen in Fig.3, where starting with small tan$`\beta `$ the cross section decreases, leading to a minimum at about tan$`\beta `$=10, and then increases again for larger tan$`\beta `$. At the minimum one has $`\sigma _{\stackrel{~}{\chi }_1^0p}1\times 10^{12}`$ when tan$`\beta `$=10 and $`m_{1/2}=600`$ GeV. More generally one has

$$\sigma _{\stackrel{~}{\chi }_1^0p}<1\times 10^{10}\mathrm{pb}$$
(6)
for the parameter domain when $`4\stackrel{<}{}\mathrm{tan}\beta \stackrel{<}{}20`$, $`m_{1/2}\stackrel{>}{}450\mathrm{G}\mathrm{e}\mathrm{V}(m_{\stackrel{~}{g}}\stackrel{>}{}1.1\mathrm{TeV})`$, $`\mu <0`$. In this domain, $`\sigma _{\stackrel{~}{\chi }_1^0p}`$ would not be accessible to any of the currently planned detectors. However, mSUGRA also then predicts that this could happen only when the gluino and squarks have masses greater than 1 TeV (and for only a restricted region of tan$`\beta `$) a result that could be verified at the LHC.

## 3 Nonuniversal SUGRA Models

In the discussion of SUGRA models with nonuniversal soft breaking, universality for the first two generations of squark and slepton masses at $`M_G`$ is usually maintained to suppress flavor changing neutral currents. One allows, however, the Higgs and third generation squark and slepton masses to become nonuniversal. We maintain gauge and gaugino mass unification at $`M_G`$.

While these models contain a large number of new parameters, their effects on $`\sigma _{\stackrel{~}{\chi }_1^0p}`$ can be charcterized approximately by the signs of the deviations from universality. One choice can greatly increase $`\sigma _{\stackrel{~}{\chi }_1^0p}`$, by a factor of 10-100 compared to the universal case, and the reverse choice can reduce $`\sigma _{\stackrel{~}{\chi }_1^0p}`$ ( though by a much lesser amount). Thus it is possible for detectors to probe regions of smaller tan$`\beta `$ with nonuniversal breaking, and detectors obeying Eq. (1) can probe part of the parameter space for tan$`\beta `$ as low as $`tan\beta 4`$.

The minimum cross section occurs (as in mSUGRA) at the lowest tan$`\beta `$ and at the largest $`m_{1/2}`$ i.e. in the coannihilation region. We limit ourselves here to the case where only the Higgs masses are nonuniversal. One finds then results similar to mSUGRA i.e. $`\sigma _{\stackrel{~}{\chi }_1^0p}\stackrel{>}{}10^9`$ pb for $`\mu >0`$, $`m_{1/2}`$ 600 GeV. For $`\mu <0`$, there can again be a cancellation of matrix elements reducing the cross section to $`10^{12}`$ pb when $`m_{1/2}=600`$ GeV in a restricted part of the parameter space when tan$`\beta `$10.

## 4 Summary

We have examined here the neutralino-proton cross section for a number of SUGRA type models. In all the models considered, there are regions of parameter space with $`\stackrel{~}{\chi }_1^0p`$ cross sections of the size that could be observed with current detectors. Thus with the sensitivity of Eq. (1), detectors would be sampling regions of the parameter space for mSUGRA where tan$`\beta \stackrel{>}{}25`$, $`m_{\stackrel{~}{\chi }_1^0}\stackrel{<}{}90\mathrm{G}\mathrm{e}\mathrm{V}`$ and $`\mathrm{\Omega }_{\stackrel{~}{\chi }_1^0}h^2\stackrel{<}{}0.1`$. Nonuniversal models can have larger cross sections and so detectors could sample down to tan$`\beta \stackrel{>}{}4`$, while for the D-brane models considered, detectors could sample down to tan$`\beta \stackrel{>}{}15`$.

The minimum cross sections these models predict are considerably below current sensitivity. Thus for mSUGRA one finds for $`\mu >0`$ that $`\sigma _{\stackrel{~}{\chi }_1^0p}\stackrel{>}{}1\times 10^9\mathrm{pb}`$ for $`m_{1/2}600\mathrm{G}\mathrm{e}\mathrm{V},\mu >0`$, where $`m_{1/2}=600`$ GeV corresponds to $`m_{\stackrel{~}{g}}1.5`$ TeV, $`m_{\stackrel{~}{\chi }_1^0}240`$ TeV. This is still in the range that would be accessible to detectors being planned (such as GENIUS or Cryoarray). For $`\mu <0`$, a cancellation can occur in certain regions of parameter space allowing the cross sections to fall below this. Thus

$`\sigma _{\stackrel{~}{\chi }_1^0p}`$ $`<`$ $`1\times 10^{10}\mathrm{pb}\mathrm{for}\mathrm{\hspace{0.17em}4}\stackrel{<}{}\mathrm{tan}\beta \stackrel{<}{}20,`$
$`\mu `$ $`<`$ $`0,m_{1/2}\stackrel{>}{}450\mathrm{GeV}`$ (7)
and reaching a minimum of $`\sigma _{\stackrel{~}{\chi }_1^0p}1\times 10^{12}`$ pb for tan$`\beta `$=10, $`m_{1/2}`$=600 GeV, $`\mu <0`$. This domain would appear not to be accessible to future planned detectors. Since $`m_{1/2}`$=450 GeV corresponds to $`m_{\stackrel{~}{g}}1.1`$ TeV, this region of parameter space would imply a gluino squark spectrum at the LHC above 1 TeV.

The above results holds for the mSUGRA model. While a full analysis of coannihilation has not been carried out for the nonuniversal and D-brane models, results similar to the above hold for these over large regions of parameter space. Thus for nonuniversal Higgs masses and for the D-brane model one finds $`\sigma _{\stackrel{~}{\chi }_1^0p}\stackrel{>}{}10^9`$ pb for $`\mu >0`$, while a cancellelation allows $`\sigma _{\stackrel{~}{\chi }_1^0p}`$ to fall to $`10^{12}`$ pb for $`\mu <0`$ at tan$`\beta 10`$ (with again a gluino/squark mass spectrum in the TeV domain).

## 5 Acknowledgement

This work was supported in part by NSF grant no. PHY-9722090.