text:
# The Nature of V359 Centauri Revealed: New Long-Period SU UMa-Type Dwarf Nova

## 1 Introduction

Cataclysmic variables (CVs) are close binary systems consisting of a white dwarf and a red dwarf secondary transferring matter via the Roche lobe overflow (for a review of CVs, see Warner (1995a)). CVs are subdivided into several categories, including dwarf novae (DNe) and novae. Both DNe and novae are characterized by the presence of a sudden increase of brightness (outburst). Although the mechanisms of DN-type outbursts (cf. Osaki (1996)) and nova outbursts (cf. Starrfield & Sparks (1987); Starrfield (1999); Starrfield et al. (2000)) are different, observational discrimination between rarely outbursting DNe and novae can be sometimes difficult (see Downes & Margon (1981) and Kato et al. (2001b) for classical and recent examples, respectively). Since rarely outbursting DNe can be easily confused with very fast novae, these confusions may have skewed our statistical view of classical novae (Downes, 1986; Liller & Mayer, 1987; Shafter, 1997).

A large fraction of such confusions turned out to be SU UMa-type dwarf novae or WZ Sge-type dwarf novae (Kato et al., 2001b). SU UMa-type dwarf novae are a subclass of DNe. WZ Sge-type dwarf novae are still enigmatic, both in theory and to observations, SU UMa-type dwarf novae, which very infrequently (once in $``$10 yr) show large-amplitude ($``$8 mag) outbursts (Bailey, 1979; Downes & Margon, 1981; Patterson et al., 1981; O’Donoghue et al., 1991). All SU UMa-type dwarf novae, including WZ Sge-type dwarf novae, show superhumps during their long, bright outbursts (superoutbursts). \[For a recent review of dwarf novae and SU UMa-type dwarf novae, see Osaki (1996) and Warner (1995b), respectively.\] Superhumps have periods (superhump period: $`P_{\mathrm{SH}}`$) a few percent longer than the orbital periods ($`P_{\mathrm{orb}}`$) (Vogt, 1980; Warner, 1985), which is believed to be a consequence of the apsidal motion (Osaki, 1985; Molnar & Kobulnicky, 1992) of a tidally induced eccentric accretion disk (Whitehurst, 1988; Hirose & Osaki, 1990; Lubow, 1991). WZ Sge-type dwarf novae are known to show a different kind of (super)humps during the earliest stage of superoutbursts (Kato et al., 1996; Matsumoto et al., 1998; Ishioka et al., 2002; Osaki & Meyer, 2002; Kato, 2002a). These (super)humps in WZ Sge-type dwarf novae have periods close to $`P_{\mathrm{orb}}`$, which can be easily distinguished from usual SU UMa-type superhumps. The presence of superhumps thus provides a powerful photometric tool in discriminating novae and SU UMa-type/WZ Sge-type dwarf novae once an object undergoes another outburst.

V359 Cen was originally discovered as a possible nova by A. Opolski (see Duerbeck (1987)). The object was visible on 19 plates taken between 1939 April 20 and 27, and the recorded maximum was mpg = 13.8 (Duerbeck, 1987). After an examination of Harvard plates of the corresponding epoch and Opolski’s finding chart, Duerbeck (1987) suggested a 21.0 mag quiescent counterpart. The true nature of the object, however, remained uncertain. The object was even proposed to be a nova in the Galactic halo. From distant nova candidates, Kato et al. (2001b) selected V359 Cen as a candidate for a rarely outbursting dwarf nova. A finding chart of the proposed quiescent counterpart was presented in Duerbeck (1987).

Munari & Zwitter (1998) tried to study the proposed quiescent counterpart spectroscopically, but the attempt failed because of its faintness ($`V`$ fainter than 20.5). Gill & O’Brien (1998) obtained a deep image around V359 Cen, and showed that the profile is indistinguishable from that of a normal star; there was no evidence of a nova shell.

The situation dramatically changed when one of the authors (Rod Stubbings) detected the second historical outburst on 1999 July 13 (vsnet-alert 3216).111 http://www.kusastro.kyoto-u.ac.jp/vsnet/Mail/alert3000/
msg00216.html.  The object further underwent outbursts in 2000 May, 2001 April and 2002 June. We photometrically observed two outbursts (1999 July and 2002 June) and revealed that V359 Cen is an SU UMa-type dwarf nova. Woudt & Warner (2001) obtained time-resolved CCD photometry following the 1999 July outburst and detected a periodicity of 0.0779 d (112 min), but interpretation of this period remained rather uncertain.

## 2 Observations

The 1999 observation by the MOA team was performed using a 61 cm Ritchey-Chrétien Cassegrain telescope (f/6.25) with the MOA-cam2 (Yanagisawa et al., 2000), constructed with three SITe back-illuminated CCDs (2047$`\times `$4095 pixels). The MOA blue filter (MOA B) covers 395–620 nm and MOA red filter covers 620–1050 nm. The exposure times were 300 and 180 s for the 1999 July 14 and 15 data, respectively. The magnitudes of the object were measured with Dophot package. The absolute calibration of the magnitudes was done using an ensemble of $``$40 neighboring stars, whose zero-point was determined using about 100 LMC standard stars measureed with the Hubble Space Telescope (HST). The MOA magnitudes can be linked to the standard $`V`$ and $`R_\mathrm{c}`$ systems using Eq. 1, where red and blue denote MOA red and MOA blue magnitudes (Noda et al., 2002). Since the blue and red observations were not completely simultaneous, we list the magnitudes on the MOA photometric system in Table 1.

$`V`$ $`=`$ $`blue0.16(bluered)+const_1`$
$`R_\mathrm{c}`$ $`=`$ $`red+0.29(bluered)+const_2`$ (1)
The 2002 observations were undertaken by the VSNET Collaboration.222 http://www.kusastro.kyoto-u.ac.jp/vsnet/. The equipment and reduction software are summarized in Table 2. The Kyoto observations were analyzed using the JavaTM-based PSF photometry package developed by one of the authors (TK). The other observers performed aperture photometry. The magnitudes were given relative to GSC 7750.220, whose constancy during the observation was confirmed by a comparison with USNO-A1.0 0450.13739601. All systems are close to $`R_\mathrm{c}`$. The journal of the 2002 observations are summarized in Table 3.

Barycentric corrections to the observed times were applied before the following analysis.

## 3 Astrometry and Quiescent Counterpart

Astrometry of the outbursting V359 Cen was performed on CCD images taken by R. Santallo (2002 June 1). An average of measurements of five images (GSC$``$2.2 system, about 20 reference stars; internal dispersion of the measurements was $``$0′′.05) has yielded a position of 11h 58m 15s.330, $``$41 46 08′′.44 (J2000.0). The position agrees with the GSC$``$2.2 star at 11h 58m 15s.322, $``$41 46 08′′.35 (epoch 1995.392 and magnitudes $`r`$ = 18.46, $`b`$ = 19.15) and the USNO$``$A2.0 star at 11h 58m 15s.330, $``$41 46 09′′.16 (epoch 1982.262 and magnitudes $`r`$ = 17.7, $`b`$ = 18.7).

This identification confirms the quiescent magnitude ($`V`$ = 18.7) reported by Woudt & Warner (2001). The quiescent magnitudes (21 or $`V`$ fainter than 20.5) reported by Duerbeck (1987) and Munari & Zwitter (1998) seem to be underestimated. An examination of POSS I red plate (limiting magnitude $`R`$ 18.5) shows that the object was near the detection limit (presumably a result of a large air-mass). This impression may have affected the estimate by Duerbeck (1987).

The failure by Munari & Zwitter (1998) in obtaining a quiescent spectrum is, however, difficult to reconcile with the value of $`V`$ = 18.7. Since a few dwarf novae categorized to established or suspected SU UMa-type dwarf novae are known to show high and low states in quiescence (HT Cas: Zhang et al. (1986); Wood et al. (1995); Robertson & Honeycutt (1996); IR Com: Richter & Greiner (1995); Kato et al. (2002a) and less established BZ UMa: Kaluzny (1986)). Although it is still premature to draw a firm conclusion, V359 Cen may belong to a small class of SU UMa-type dwarf novae with high/low transitions in quiescence.

## 4 Long-Term Light Curve

Fig. 1 shows the long-term visual light curve of V359 Cen constructed from the observations reported to the VSNET Collaboration. Large and small dot represent positive and negative (upper limit) observations, respectively. Four outbursts (1999 July, starting on JD 2451373; 2000 May, on JD 2451680; 2001 April, on JD 2452025 and 2002 June, on JD 2452422) are clearly seen. The intervals between the detected outbursts are in the range of 307–397 d. Although there were unavoidable seasonal gaps in observations, these values seem to be a representative outburst cycle length. The observed maximum magnitudes of the outbursts were $``$13.8. This constancy of the maximum magnitudes likely precludes the previously supposed possibility that the maximum of the 1939 outburst was missed (Duerbeck, 1987).

## 5 The 1999 Outburst

Fig. 2 shows the enlarged light curve of V359 Cen during the 1999 outburst. The data are from the MOA observations. Although complete phase coverage was impossible because of the short available runs, the object clearly exhibited variations with amplitudes of $``$0.15 mag, which can be attributed to superhumps.

A period analysis of the MOA blue data using Phase Dispersion Minimization (PDM; Stellingwerf (1978)) has yielded the theta diagram presented in Fig. 3. Although a unique alias selection is impossible from these data only, we can safely choose the correct alias of $`P_{\mathrm{SH}}`$ = 0.0824(4) d based on the later determination of the superhump period (Sect. 6.2).

## 6 The 2002 Outburst

### 6.1 Course of the Outburst

The 2002 outburst was detected by Rod Stubbings on May 28 at $`m_{\mathrm{vis}}`$ = 13.8 (vsnet-alert 7356).333 http://www.kusastro.kyoto-u.ac.jp/vsnet/Mail/alert7000/
msg00356.html.  CCD time-resolved photometry started within a day following this detection. Fig. 4 shows the overall light curve of the 2002 superoutburst. The magnitudes are relative to GSC 7750.220 and are on a system close to $`R_\mathrm{c}`$. After two days of the outburst, the system started to fade linearly at a rate of 0.16 mag d-1. This slowly and linearly fading phase (often referred to as superoutburst plateau) is very characteristic of an SU UMa-type superoutburst (Vogt, 1980; Warner, 1985).

The mean decline rate of 0.16 mag d-1 during the plateau phase is larger than those of other SU UMa-type dwarf novae (Table 4) with similar $`P_{\mathrm{SH}}`$ to that of V359 Cen (cf. Sect. 6.2). Since a higher mass-transfer rate from the secondary star tends to thermally stabilizes the accretion disk and reduce the decay rate (e.g. Osaki (1995)), a rather exceptionally large decay rate in V359 Cen may be a result of a systematically smaller mass-transfer rate.

### 6.2 Superhump Period and Evolution

Fig. 5 shows enlarged nightly light curves. Superhumps are clearly visible on all observed nights. Fig. 6 shows the result of a PDM period analysis of the 2002 data between May 31 and June 6 (superoutburst plateau). The linear declining trend has been subtracted before the analysis. The resultant best $`P_{\mathrm{SH}}`$ is 0.08092(1) d. The selection of the correct alias has been confirmed by independent period analyses of individual long continuous runs. This period established that V359 Cen is a long-period SU UMa-type dwarf nova. Fig. 7 shows the phase-averaged profile of superhumps. The rapidly rising and slowly fading superhump profile is characteristic to an SU UMa-type dwarf nova (Vogt, 1980; Warner, 1985).

We extracted the maxima times of superhumps from the light curve by eye. The averaged times of a few to several points close to the maximum were used as representatives of the maxima times. The errors of the maxima times are usually less than $``$0.002 d. The resultant superhump maxima are given in Table 5. The values are given to 0.0001 d in order to avoid the loss of significant digits in a later analysis. The cycle count ($`E`$) is defined as the cycle number since BJD 2452423.939. A linear regression to the observed superhump times gives the following ephemeris:

$$\mathrm{BJD}(\mathrm{maximum})=2452423.9503+0.08108E.$$
(2)
Fig. 8 shows the ($`OC`$)’s against the mean superhump period (0.08108 d) from Eq. 2. Although the general trend can be expressed by a negative quadratic term of $`\dot{P}`$ = $`10.8(1.5)\times 10^6`$ d cycle-1, or $`\dot{P}/P`$ = $`13.3(1.9)\times 10^5`$, superhump maxima with $`E23`$ can be well expressed by a constant period of $`P`$ = 0.08094 d and $`|OC|`$’s less than 0.005 d. This finding indicates that the superhump period was virtually constant during the plateau phase (the nominal $`\dot{P}`$ during this period os $`5.3(2.8)\times 10^6`$ d cycle-1). A sudden change between $`E`$ = 0 and $`E`$ = 23 can be interpreted as a result of rapid evolution of superhumps during the earliest stage of the superoutburst. The “textbook” evolutionary time-scales (2–3 d) of superhumps in long-period SU UMa-type dwarf novae (Warner, 1985) also support this interpretation. The large period change of $`P_{\mathrm{SH}}`$ observed during the earliest stage of the 2002 superoutburst may explain the slight discrepancy of the periods between the 1999 and 2002 observations (the 1999 observation corresponds to the earliest stage of a superoutburst). The period determined from $`0E23`$ observations is 0.0817 d, sufficiently longer than the most likely orbital period (Sect. 7.2), precludes the possibility of WZ Sge-type early (super)humps as the origin of these early modulations.

### 6.3 Super-QPOs

Some SU UMa-type dwarf novae are known to show large-amplitude, highly coherent quasi-periodic oscillations (QPOs) during the evolution stage of superhumps (Kato et al., 1992; Kato, 2002b). These QPOs are sometimes referred to as “super-QPOs”. The light curve of V359 Cen on May 29 shows a hint of such QPOs. The lower panel of Fig. 9 shows residuals of the May 29 light curve after subtracting the superhump signal by using a Fourier decomposition of the superhump profile up to the third harmonics, and subtracting a slow linear trend. Small-amplitude modulations with a typical time scale of $``$0.02 d are present. Fig. 10 shows a power spectrum of the residual shown in the lower panel of Fig. 9. The strongest power is present at a frequency of $``$49 d-1, which corresponds to a period of $``$0.02 d.

Although the QPOs were not as prominent as seen in SW UMa (Kato et al., 1992) or (Kato, 2002b), the absence of similar signals on later nights suggest that this variation is a kind of super-QPOs. As shown in Sect. 6.2, the epoch of the detection of QPOs corresponds to the rapidly evolving stage of superhumps. This finding further supports an idea that super-QPOs are associated with the growth of superhumps (Kato et al., 1992; Kato, 2002b)

## 7 V359 Cen as an SU UMa-Type Dwarf Nova

### 7.1 Outburst Characteristics

Table 6 lists the observed parameters of the four recorded outbursts. All the recorded outbursts have durations longer than 4 d, which indicate that all the recorded outbursts were superoutbursts. The intervals between these outbursts (307–397 d) are thus regarded as the supercycle of this SU UMa-type dwarf nova. Since the present duty cycle of observations is about 50%, there remains a small possibility that the true supercycle is the half of this value. Further dense monitoring is strongly encouraged to completely exclude this possibility. It is well established that almost all well-studied SU UMa-type dwarf novae with similar supercycle lengths show a few to $``$ten normal outburst within one supercycle (Warner, 1995b; Nogami et al., 1997). Since normal outbursts of SU UMa-type dwarf novae are usually only 0.5–1.0 mag fainter than superoutbursts (Warner, 1985), some of such normal outbursts should have been detected. The present non-detection of normal outbursts seems to suggest that V359 Cen has fewer normal outbursts than in usual SU UMa-type dwarf novae. Although future deeper monitoring for normal outbursts is absolutely needed, the occurrence of normal outbursts in V359 Cen may be effectively suppressed by an unknown mechanism; similar instances have been reported in other SU UMa-type dwarf novae (Kato, 2001; Kato et al., 2002b). If normal outbursts are entirely missing, V359 Cen would become the first long-period analog of V844 Her which has been known to show only superoutbursts (Kato & Uemura (2000): $`P_{\mathrm{SH}}`$ = 0.05592 d, supercycle lengths = 220–290 d). Similar SU UMa-type dwarf novae with predominance of superoutbursts are highly concentrated in the short period region (cf. Kato & Uemura (2000)). V359 Cen may be an exceptional object in its combination of outburst characteristics and the superhump (or orbital) period.

We also note that the outburst characteristics of V359 Cen also resemble those of long-period SU UMa-type dwarf novae EF Peg (Howell et al., 1993; Kato, 2002b) and V725 Aql (Uemura et al., 2001). Both EF Peg and V725 Aql only infrequently show normal outbursts, which is exceptional among long-$`P_{\mathrm{orb}}`$ SU UMa-type dwarf novae (cf. Warner (1995b)). While long outburst recurrence times imply that these systems have low mass-transfer rates (Ichikawa & Osaki, 1994; Osaki, 1996), recent detailed calculations of the evolution of CVs (e.g. Podsiadlowski et al. (2001)) suggest that mass-transfer rates are higher in long-period systems even if the effect of stellar core evolution is properly taken into account. These systems (EF Peg, V725 Aql and possibly V359 Cen) may be violating the modern evolutionary scenario of CVs. Future determination of the binary parameters and stellar composition analysis of these systems are therefore strongly encouraged.

### 7.2 Superhump Excess and Late Superhumps

Woudt & Warner (2001) tentatively identified their photometric period (112 min) as late superhumps (Haefner et al., 1979; Vogt, 1983; van der Woerd et al., 1988; Hessman et al., 1992), which are known to have similar periods with ordinary superhumps (i.e. a few percent longer than $`P_{\mathrm{orb}}`$), but have phases of $``$ 0.5 different from those of ordinary superhumps. Since this signal was observed long after then cessation of the 1999 outburst, Woudt & Warner (2001) suggested that V359 Cen may have shown a long persistence of late superhumps as was observed in a WZ Sge-type dwarf nova in EG Cnc (Kato et al., 1997; Patterson et al., 1998).

The present correct identification of the superhump period, however, indicates that the 112 min periodicity observed by Woudt & Warner (2001) can not be attributed to superhumps, but can be better understood to represent $`P_{\mathrm{orb}}`$. This observation, on the contrary to the suggestion by Woudt & Warner (2001), indicates that the superhumps or late superhumps in V359 Cen must have decayed more rapidly, as in other usual SU UMa-type dwarf novae (see Kato et al. (2001a) for a recent example of the decay of late superhumps). V359 Cen is thus unlikely related to WZ Sge-type dwarf novae which always show long persistence of late superhumps.

By adopting $`P_{\mathrm{orb}}`$ = 0.0779 d, we obtain a fractional superhump excess ($`ϵ=P_{\mathrm{SH}}/P_{\mathrm{orb}}1`$) of 3.9% for the best $`P_{\mathrm{SH}}`$ is 0.08092(1) d (cf. Sect. 6.2). Using the “mean” $`P_{\mathrm{SH}}`$ = 0.08108 d from the entire 2002 superoutburst, we obtain $`ϵ`$ = 4.1%. These fractional superhump excesses are not unusual for an SU UMa-type dwarf nova with this $`P_{\mathrm{orb}}`$ (Molnar & Kobulnicky, 1992; Patterson, 1998), suggesting that V359 Cen should have a normal binary mass ratio $`q=M_2/M_1`$ in spite of its rather unusual outburst characteristics.

## 8 Summary

We detected four outbursts of V359 Cen (possible nova discovered in 1939) between 1999 and 2002. Time-resolved CCD photometry during two outbursts (1999 and 2002) revealed that V359 Cen is actually a long-period SU UMa-type dwarf nova with a mean superhump period of 0.08092(1) d. We identified its supercycle length as 307–397 d. This secure identification of the superhump period precludes the previously supposed possibility that V359 Cen could be related to a WZ Sge-type system with a long persistence of late superhumps. The outburst characteristics of V359 Cen are, however, rather unusual in its low occurrence of normal outbursts. The fractional superhump excess is 3.9–4.1%, which suggests that V359 Cen should have a normal binary mass ratio in spite of its rather unusual outburst characteristics. We also obtained a secure identification of the quiescent counterpart and discussed on the possibility of high/low state changes. The evolution of superhumps and their period change was closely followed. We also detected super-QPO-type variation (period $``$0.02 d) during the earliest stage of the 2002 superoutburst.

The authors are grateful to Andrew Pearce who reported visual observations of V359 Cen to VSNET. This work is partly supported by a grant-in aid (13640239) from the Japanese Ministry of Education, Culture, Sports, Science and Technology. Part of this work is supported by a Research Fellowship of the Japan Society for the Promotion of Science for Young Scientists (MU). The CCD operation of the Bronberg Observatory is partly sponsored by the Center for Backyard Astrophysics. The CCD operation by Peter Nelson is on loan from the AAVSO, funded by the Curry Foundation. This research has made use of the Digitized Sky Survey producted by STScI, the ESO Skycat tool, the VizieR catalogue access tool, and the USNOFS Image and Catalogue Archive operated by the United States Naval Observatory, Flagstaff Station (http://www.nofs.navy.mil/data/fchpix/).


format: text

text:
# Dynamical mean-field approximation for pair contact process with a particle source

The study of phase transitions from active fluctuating phases into absorbing states has attracted considerable interest in the last decade . The static exponents of these phase transitions belong frequently to the directed percolation (DP) universality class, however, there are also examples for new universality classes . Open questions are related to the necessary conditions that can destroy this robust universal behavior. One of the simplest model showing static DP behavior is the pair contact process (PCP) which was proposed to realize a system involving infinitely many absorbing states . Recently, Dickman et. al. introduced a modified PCP model to explore the robustness of DP transition . In this model, each site of the one-dimensional lattice is either vacant or occupied by a single particle. A (randomly chosen) pair of nearest-neighbor particles is annihilated with a probability $`p`$ or an additional particle is created around the given pair with a probability $`1p`$ if it is not forbidden by double occupancy. In the extended model, an external particle source is introduced that attempts to insert isolated particles with a rate of $`h`$. This system exhibits an active phase when $`p`$ is smaller than a critical value, $`p_c`$. For $`pp_c`$ the system evolves into a frozen (absorbing) state where the nearest-neighbor pairs are absent. Further details of the model can be found in Ref. . This extended model was studied by Monte Carlo (MC) simulations and dynamical cluster mean-field approximation for quite large cluster sizes (ranged from $`n=2`$ to $`6`$). Some disturbing behaviors, however, remained unsolved. For example, the analytical predictions tend not monotonously toward the MC results for $`h>0`$ when the cluster size is increased. Furthermore the analytical results shows a discontinuity in the variation of critical point $`p_c`$ if $`h0`$. This observation is surprising because the present approximation has proved to be satisfactory in many cases for $`n6`$ . It is expected that the further increase of cluster size will resolve these discrepancies.

In this Brief Report, we discuss the results of dynamical mean-field approximations for cluster sizes as large as $`n=12`$. The present approach allows us to give more accurate extrapolations to the limit $`n\mathrm{}`$. Using this large number of data (from $`n=212`$) we can also improve the prediction of coherent anomaly method (CAM)  when considering the critical exponents around the transition point. Furthermore, we could study the frozen density of isolated particle concentration toward the system evolves in the inactive phase .

To identify the transition point we have evaluated the $`p`$-dependence of the density of pairs of particles $`\rho _2`$ (henceforth considered as order parameter) as well as the density of isolated particles $`\rho _1`$. The transition point is determined by the zero-point of $`\rho _2`$ for $`h=0`$ and by the breaking point of the function $`\rho _1(p)`$ for $`h>0`$. This latter criteria was particularly useful because of the extremely low value of $`\rho _2`$ in a large interval of $`p`$.

We shall not describe the details of the mean-field technique because it has already been applied and demonstrated for the PCP model by several authors previously . The derivation of the hierarchy of equations of motion for the configuration probabilities on the $`n`$-site clusters becomes complicated for large $`n`$. The technical difficulties increase drastically when enlarging the size of clusters. Traditionally, the set of equation of motion is solved numerically in the stationary states. Despite the slow convergency toward the stationary solution for large $`n`$, the numerical integration of the master equations seems to be a more efficient method to find the stationary probability of configurations than the traditional Newton-Raphson method. Using the numerical integration method we could determine all the configuration probabilities appearing on the 12-site cluster. (As an example, to get a data in last row of Table I requires 4-week running on a personal computer.) It is worth mentioning that very recently this method is used successfully for the consideration of a stochastic sandpile model .

Our results are plotted in Figs. 13. Table I summarizes the predictions of all levels of approximations for $`p_c`$ and for the density of isolated particles $`\varphi _{nat}`$ at the critical point at different values of $`h`$. Figure 1 demonstrates clearly that the predictions for $`p_c`$ become monotonously convergent for $`n6`$. The linear fit for $`n>5`$ data gives $`p_c=0.075`$ in close agreement with the MC data ($`p_c^{MC}=0.077`$). The same good convergence can be obtained in the presence of external source. The comparisons of linear fit with MC data are also listed in Table I (MC data are taken from Ref ).

Surprisingly, although the $`n`$-site approximations converge to MC data at arbitrary value of $`h`$, a jump can be observed in $`p_c`$ when considering the limit $`h0`$ for all $`n`$. The magnitude of this jump is defined as

$$\mathrm{\Delta }p^n=\underset{h0}{lim}p_c^n(h)p_c^n(h=0).$$
(1)
The inset of Fig 1 shows $`\mathrm{\Delta }p^n`$ for different levels of the approximation. This log-log plot suggests that the jump decreases as a power low function, i.e., $`\mathrm{\Delta }p^nn^\omega `$ with an exponent $`\omega 1.666`$. The explanation of this value of $`\omega `$ remains to be clarified. This discontinuity of $`\mathrm{\Delta }p^n`$ may also be observed for higher dimension versions of the model .

A distinct improvement of convergency can also be observed for higher values of $`n`$ if we consider the density of isolated particles $`\varphi _{nat}`$ at $`p_c`$. The extrapolations, based on predictions of $`n>5`$ data, are in excellent agreement with MC data.

Such a large number of approximations makes possible to apply the CAM analysis introduced by Suzuki . In the vicinity of the critical point the order parameter ($`\rho _2`$) and the density of isolated particles ($`\rho _1`$) are estimated by $`\rho _2a_2^n(p_c^np)`$ and $`\rho _1\varphi _{nat}a_1^n(p_c^np)`$, where $`p_c^n`$ denotes the prediction for the critical point at $`n`$-point level. To estimate the $`\beta `$ critical exponents we have plotted the amplitudes of mean-field results ($`a_1^n`$ and $`a_2^n`$) as a function of

$$\delta _n=(p_c^n/p_c)^{1/2}(p_c/p_c^n)^{1/2},$$
(2)
where $`p_c`$ denotes the result of MC simulation. Figure 2 shows the mean-field amplitudes in the absence of external source ($`h=0`$) for different levels of approximation.

Using different subsets of CAM-data the estimated exponents are $`0.221,0.246,0.276,`$ and $`0.312`$. The dashed lines correspond to the $`\beta =0.276`$ DP exponent. Although the error bar is still large but the tendency to DP is straightforward. These results confirm the earlier MC observations . Namely, the order parameter and the density of isolated particles can be described by the same power laws when the system approaches the critical point from the active phase. At the same time the CAM analysis of mean-field data helps to understand why the exponents agree. Figure 2 demonstrates that these two quantities are proportional to each other in the vicinity of critical point. (In this plot $`1.1a_1^n`$ is compared with $`a_2`$).

In the presence of source ($`h>0`$) the earlier MC simulation suggested a slightly modified $`\beta =0.287`$ exponent . Our mean-field data become also convergent for $`n6`$ suggesting an exponent close to DP class ($`\beta =0.274\pm 0.038`$ at $`h=10`$). However, the region in $`\delta _n`$ is so small that the error bar of the estimated exponent is four times larger than the predicted change of MC result. Therefore, this method is incapable to distinguish such close exponents.

As noted above, the large-cluster approximation makes possible to test the prediction of static critical behavior of particle concentration in the inactive phase . Such a critical behavior of a static quantity in the inactive phase has already been observed in a different model . In a previous MC study of PCP model it was suggested that the density of isolated particles follows a power law $`\rho _1^{nat}\rho _1(pp_c)^{\beta _1}`$ for $`p>p_c`$ in the absence of external source . Following Marques et al. , we have calculated the absorbing state selected by the system’s dynamics starting from the same homogeneous initial state. However, different initial particle concentrations were chosen to study the robustness of the final absorbing state in the inactive phase. The initial concentration is ranged from $`\varphi _n`$ to $`0.99`$. The results of $`6`$-point approximation are plotted in Fig. 3 for different values of initial concentration at $`h=0`$. It suggests that the absorbing state selected by the system’s dynamics depends on the initial condition. Similar behavior may be observed for all levels of approximation. The observed deviation of curves increases further by increasing $`p`$ and the level of approximation. Obviously, the stationary value of concentration becomes independent of the initial condition in the active phase.

We have tried to apply CAM analysis on the data obtained from the same initial concentration for different levels of approximation. Even in this case, the irregularity of mean-field coefficients at different levels of approximation does not allow to extract critical exponent. Shortly, the present approximation does not support the existence of static critical behavior in the inactive phase. In the light of this prediction further intensive MC simulations are suggested to clarify the possible existence of natural absorbing states.

In summary, we have demonstrated that the present dynamical mean-field approximations (for large $`n`$) yield adequate extrapolation to the limit $`n\mathrm{}`$ even in the presence of external particle source that causes nonanalytical behavior in the limit $`h0`$. According to this approximation the discontinuity of $`p_c`$ displays power-law decay. The application of CAM analysis demonstrates that the densities of nearest-neighbor pairs and isolated particles are proportional to each other in the vicinity of the critical point and supports the same critical exponents of these quantities. This behavior is certainly related to the fact that the extinct pairs can leave extra solitary particles behind. Although our analysis is proved to be a useful tool for investigating stationary states in the active phase it does not confirm the occurrence of an unambiguous critical behavior in the inactive phase where the composition of the frozen state depends on the initial state. We hope that the obscurity of stationary inactive state will stimulate further MC simulations and theoretical investigations.

The author wishes to thank Ron Dickman, György Szabó, and M. C. Marques for stimulating discussions. This research was supported by the Hungarian National Research Fund (OTKA) under Grant No. F-30449 and Bolyai(BO/0067/00).

Electronic address: szolnoki@mfa.kfki.hu


text:
# NEUTRINO PHYSICS AND ASTRONOMY WITH MACRO

## 1 MACRO as an atmospheric neutrino detector

The MACRO detector was located in the Gran Sasso underground laboratory . Thanks to its large area, fine tracking granularity and up-down symmetry, it was a proper tool for the study of upward-travelling muons and neutrino interactions in the apparatus. The different kinds of neutrino events detected by MACRO are shown in Fig. 1A : (1) upward-throughgoing muons, (2) semicontained upgoing muons, (3) upgoing stopping and (4) semicontained downgoing tracks. The sample (1) is due to more energetic neutrinos ($`<E_\nu >50GeV`$) producing muons also at long distances from the detector. The other samples are due to neutrinos of lower energy ($`<E_\nu >4GeV`$). The samples (3) and (4) are indistinguishable and therefore they are studied together (3+4).

In Fig.s 1B and 2 the angular measurements are shown for samples (1), (2) and (3+4), respectively. All distributions are so far from the expectation assuming no oscillation. Assuming $`\nu _\mu \nu _\tau `$ oscillation the measurements become compatible with expectation . The best fit parameters result $`\mathrm{sin}^22\theta _{mix}=1`$ and $`\mathrm{\Delta }m^2=0.0025eV^2`$. The measurement of the Multiple Coulomb Scattering  permits to estimate the energy of muons in sample (1). Neutrino energies are inferred by means of Monte Carlo methods. In Fig. 3A the data/expectation ratio as a function of estimated $`L/E_\nu `$ is shown. The last point is due to sample (2).

In Fig. 3B the 90% C.L. allowed regions assuming $`\nu _\mu \nu _\tau `$ oscillation are shown. The smaller area is estimated by normalization and angular distribution of sample (1). The medium area is due to $`\mu `$-energy estimate for the same sample. The larger area is deduced by (2) and (3+4) low energy samples.

## 2 Neutrino astronomy

A sample of 1356 upward-travelling muons ($`E_\mu >1GeV`$) has been studied to look for astrophysical signals . The data do not show $`\nu `$-excess from selected pointlike sources. Neutrino-induced muon-flux limits have been established at the level of $`10^{14}cm^2s^1`$. Recently the microquasar GX339-4 has been proposed as a source of TeV neutrinos . We note that MACRO detected 7 events in the direction of GX339-4, while 2.3 atmospheric neutrinos were expected . Another analysis searched for high-energy upward-travelling muons due to astrophysical high-energy diffuse neutrinos. No statistically significant signal has been found and an upper limit has been set on diffuse neutrino flux assuming $`E_\nu ^2`$ as power law spectrum .


text:
# Canonical aspects of strangeness enhancement

## 1 INTRODUCTION

The enhanced production of strange particles in nucleus-nucleus (AA) relatively to proton-proton (pp) or proton-nucleus (pA) collisions, was since long, argued  to be a signal of a Quark-Gluon Plasma (QGP) formation. This enhancement is seen in all experiments at all energies from AGS to SPS up to RHIC. Furthermore, enhancement of strange baryons and anti-baryons was predicted to depend on their strangeness content and to appear in a typical hierarchy, $`E_\mathrm{\Lambda }<E_\mathrm{\Xi }<E_\mathrm{\Omega }`$, which was observed in particular by the WA97 , NA57  and NA49  Collaborations. These enhancement features and hierarchy, are argued to be expected at all energies as a consequence of a canonical thermal phase-space suppression.

## 2 CANONICAL STRANGENESS SUPPRESSION

When describing a statistical system there are at least two ways of implementing the conservation laws of quantum numbers. In the Grand Canonical (GC) formulation conservation of a quantum number is ensured on average by use of the corresponding chemical potential. This description is only valid when particles, with a given conserved charge, are produced with a large multiplicity (system of a large volume and/or high temperature). In the Canonical formulation (C) conservation of a quantum number is implemented exactly on an event-by-event basis. This approach is relevant when particles carrying a given conserved charge are produced with a small multiplicity (system of small volume and/or low temperature). In the usual definition of SE one compares the yield per participant of a given particle $`i`$ with strangeness $`s`$ in the large AA system to the yield per participant, of the same particle, in the small pp or pA system. Assuming that the volume parameter scales with the number of participants, $`VN_{part}`$, the enhancement of this particle is given by

$$E_s^i=\frac{\left(n_s^i\right)^{AA}}{\left(n_s^i\right)^{pp}}$$
(1)
where $`n_s^i`$ is the number density of particle $`i`$. From the above considerations, the pp (or pA) system has to be described in the canonical formulation. In principle, all conserved charges (baryon number $`B`$, electric charge $`Q`$, and strangeness $`S`$) have to be treated canonically. However, with a good approximation, $`B`$ and $`Q`$ can be treated grand canonically even in high energy pp collisions. This results in a slight overestimate of the canonical effect .

The canonical partition function of a gas of particles and resonances having strangeness $`s=0,\pm 1,\pm 2,\pm 3`$, and with total strangeness $`S=0`$ can be written as

$$Z_{S=0}^C(T,V)=\frac{1}{2\pi }_\pi ^\pi 𝑑\varphi \mathrm{exp}\left(\underset{n=3}{\overset{3}{}}S_ne^{in\varphi }\right)$$
(2)
where $`S_n=V_kZ_k^1`$ , is the sum over all particles and resonances carrying strangeness $`n`$ and where $`Z_k^1(g_k/2\pi ^2)m_k^2TK_2(m_k/T)\mathrm{exp}(B_k\mu _B/T+Q_k\mu _Q/T)`$ is the one-particle partition function of particle $`k`$ with multiplicity $`g_k`$, mass $`m_k`$, baryon number $`B_k`$, electric charge $`Q_k`$ and the corresponding chemical potentials $`\mu _B`$ and $`\mu _Q`$. From Eq.(2) one can derive  the canonical density $`(n_s^i)^C`$ and show that

$$(n_s^i)^C=(n_s^i)^{GC}F_s(T,V)$$
(3)
where $`(n_s^i)^{GC}`$ is the grand canonical density and $`F_s(T,V)`$ is the suppression factor that measures a deviation of $`(n_s^i)^C`$ from its asymptotic, grand canonical value. For large $`V`$ and/or $`T`$, $`F_s(T,V)1`$, and $`(n_s^i)^C`$ reaches its grand canonical value.

From Eq.(1), SE for a given $`VN_{part}`$ relative to pp is determined from

$$E_s^i=\frac{F_s(T,V)}{F_s(T,V_{pp})}\frac{F_s(T,N_{part})}{F_s(T,2)}$$
(4)
Figure 1 displays an example of the suppression factor at $`\sqrt{s}=8.73`$ GeV. The enhancement hierarchy at corresponding energy, as seen in Figure 2, is a direct consequence of the behavior of the suppression factor with the strangeness content of the particle . Thus, SE can indeed appear as a canonical suppression effect.

## 3 ENERGY DEPENDENCE OF STRANGENESS ENHANCEMENT

We have studied  the energy dependence of SE at four energies:$`\sqrt{s}=8.73,12.3,17.3`$ GeV (NA49, NA57, WA97, SPS) and $`\sqrt{s}=130`$ GeV (RHIC). At 17.3 and 130 GeV the values of $`T`$ and $`\mu _B`$ were taken from the thermal fit from . At $`\sqrt{s}=`$ 8.73 and 12.3 GeV the values of these parameters were determined following unified freeze-out curve  and by using the measured $`<\pi >/<N_{part}>`$ ratio at the relevant energies. Moreover, as a first approximation, the assumption have been made that the thermal parameters are independent of centrality. We have shown that SE is a decreasing function of collision energy . The result for $`\mathrm{\Lambda }`$ is given in Figure 3. In Figure 4 the model predictions are compared to the recent experimental data of NA49 Collaboration . The experimental points are deduced from the $`(<\mathrm{\Lambda }>/<\pi >)^{AA}/(<\mathrm{\Lambda }>/<\pi >)^{pp}`$ ratio  multiplied by $`<\pi >/<N_{part}>`$ in AA and in pp collisions  at the same energy. The solid line in Figure 4 is the prediction of the canonical statistical model.

## 4 CONCLUDING REMARKS

Canonical suppression plays an important role in SE. The present canonical statistical model naturally describes  the bulk properties of the WA97 data for $`\mathrm{\Lambda },\mathrm{\Xi },\mathrm{\Omega }`$, that show SE hierarchy, and saturation of the yields per participant for large $`N_{part}`$. However, the quantitative understanding of the recent data of the NA57 Collaboration , showing a smooth increase of that yields with $`N_{part}>100`$ and a significant change in the behavior of $`\overline{\mathrm{\Xi }}`$ ( the yield rises from $`N_{part}=62`$ to $`N_{part}=121`$ by a factor 2.6) would required a further study on the behavior of thermal parameters with centrality . These new features of NA57 data may be understood admitting a non-linear dependence of the canonical correlation volume or of the off-equilibrium strangeness fugacity with centrality . Nevertheless, the enhancement hierarchy and the decrease of enhancement with increasing collision energy are a generic features of the canonical statistical thermal model that are independent on the particular choice of the parameters. The most recent NA49 data on $`\mathrm{\Lambda }`$ yields seem to support the model prediction that SE is a decreasing function of collisions energy. However, experimental confirmation of these behavior by even more spectacular variation of $`\mathrm{\Xi }`$ and $`\mathrm{\Omega }`$ enhancement with $`\sqrt{s}`$ would be required. We note that UrQMD  and DPM  predicts a larger SE at RHIC than at SPS energy, which is in contrast to present statistical model results.

Acknowledgement: We wish to thank P. Braun-Munzinger, J. Cleymans, M. Gazdzicki, H. Oeschler, H. Satz and R. Stock for interesting discussions. K.R acknowledges the support of the Alexander von Humboldt Foundation.


text:
# Reply to Comment on “ Universal Fluctuations in Correlated Systems”

Bramwell et al. reply

Watkins Chapman and Rowlands claim to show that the deviation from the FTG distribution, due to the correlations between random variables introduced in the extremal statistics problem in  is the result of the slow convergence of the PDF for extreme values, with system size, for effective Gaussian variables. We think this result is correct, but it must be put in the context of recent developments in the field.

Our original motivation for studying correlated extremal statistics explicitly excludes the kind of slow relaxation towards an asymptotic, or thermodynamic limit function discussed in  and : The PDF for order parameter fluctuations in the low temperature phase of the 2D-XY model (BHP) is a thermodynamic limit function that is different to the FTG distribution. The model studied is diagonalisable into statistically independent variables, for which the dispersion in amplitudes diverges with system size. The PDF for the extreme (largest) value of these variables is the FTG distribution not the BHP distribution . Conclusion: simple extremal statistics do not explain the observed results. If extremal statistics are relevant then they must apply to more complex (correlated) many body objects, rather than the statistically independent variables of the problem. If this is not the case then extremal statistics are irrelevant for at least this model correlated system.

Our first attempt to look at the extremal statistics of correlated random variables was the model presented in  and discussed in the Comment. The authors are probably correct to conclude that the main effect comes from finite size corrections, rather than the correlations introduced in the model. This can therefore be classified as ”weak correlation” and is the extreme value equivalent of introducing a finite correlation length in a thermodynamic system and then taking the thermodynamic limit. In a strong correlation limit one would expect deviations from FTG to remain on taking the limit $`N\mathrm{}`$. This scenario has been seen in detail in reference , where the extreme value for avalanche sizes in the Sneppen depinning model is seen to follow the BHP distribution over a large range of time and length scales. We also note that renormalization group analysis of extreme value statistics for long range correlated signals shows that the tail of the resulting distribution renormalizes from the $`\mathrm{exp}(y)`$ asymptote of the FTG distribution to $`y\mathrm{exp}(y)`$ , in agreement with the exact asymptotic form for the BHP function . The results of  and  go somewhere towards confirming our hypothesis, first proposed in , that deviations to the FTG distribution introduced by correlations could provide the desired link between extremal statistics and the fluctuations of a global quantity in such correlated systems. We suggest that the results presented here (and therefore in ) and in reference  pass asymptotically close to this correlated regime.

Finally we remark that, even in the case where global fluctuations are described by the FTG distribution, as in ref  (see also ref. ), any connection with extremal statistics remains unproven. The global quantity is a sum over a macroscopic number of elements, while extremal statistics would require the selection of the biggest element by a “Maxwell’s demon”. The connection between these processes is an open problem. The authors of this comment have presented arguments based on slow relaxation towards the FTG distribution for Gaussian variables  which could be a step in the right direction. This is a useful contribution but the arguments presented here and in the references below, show that it is far from the complete picture and that it would be quite wrong to infer that deviation from the FTG distribution in the ensemble of systems discussed is generically due to finite size corrections to the thermodynamic limit. More work needs to be done.


text:
# CO Luminosity Functions For FIR and B-band Selected Galaxies and the First Estimate for Ω_{𝐻𝐼+𝐻₂}

## 1 Introduction

The star formation history of the Universe is closely linked with the evolution of the gas content of the Universe. Observations of the distribution and total gas content in galaxies and in intergalactic clouds offer some of the most important observational constraints for the cosmology and galaxy evolution models. Because stars form out of cold, dense gas clouds, the history of galaxy formation and evolution is also the history of gas accretion and conversion into stars. Here we derive the total cold gas density for the local volume by deriving the local CO luminosity function and combining with the existing estimates of the neutral atomic gas density.

A variety of luminosity functions (LFs) and mass functions are found in the literature: optical (e.g. Marzke, Huchra, & Geller, 1994), HI (e.g. Zwaan et al., 1997; Schneider et. al, 1998), and infrared (e.g. Soifer et al., 1987; Yun et al., 2001). While several groups have investigated the local neutral hydrogen density, little information is available on the molecular gas content in galaxies because of the lack of appropriate data. In her pioneer work, Verter (1987) derived averaged values of CO luminosity and the CO/HI flux ratio based on a maximum likelihood probability distribution for 40 galaxies and 47 upper limits obtained by combining all available data at the time. Since the large fraction of the CO measurements used were non-detections (upper limits), the derivation of the CO luminosity function was highly problematic. Nevertheless, a general trend of increasing number of CO emitting galaxies per unit volume with decreasing CO luminosity was suggested by this analysis.

Nearly 15 years later, the availability of the extragalactic CO data has improved greatly. The FCRAO (Five College Radio Astronomy Observatory) Extragalactic CO Survey (Young et al., 1995) represents a particularly rich database for the investigation of CO luminosity function and molecular gas content in galaxies. We have constructed a large, statistically significant sample of far-infrared (FIR) and optical $`B`$-band selected galaxies from this survey and derived a non-parametric CO luminosity function. Our sample galaxies range over 4 orders of magnitudes in CO luminosity, and only a small fraction are non-detections. Since CO is a tracer of hydrogen molecules, molecular hydrogen mass can be derived from CO luminosity. By integrating the resulting molecular gas mass function, we then derive the molecular gas mass density and the total cold gas mass density of the local universe.

## 2 Sample Selection and Properties

The FCRAO Extragalactic CO Survey (“the Survey”, herein) is the largest CO survey of galaxies available in the literature, containing observations of 300 galaxies at 1421 positions (Young et al., 1995). Since the Survey covers a broad range of Hubble types and contains a large number of detections (236 or 79%), it is well suited for a statistical analysis and for constructing a CO luminosity function. The majority of the survey galaxies were selected from the Second Reference Catalog (RC2; de Vaucouleurs et. al, 1976) or the IRAS database. Most of the galaxies in the Survey are spiral or irregular galaxies at declination north of $`25^{}`$ and satisfy at least one of the following three criteria: $`B_T^0<13`$, $`S_{60\mu m}>5`$ Jy, or $`S_{100\mu m}>10`$ Jy (see Young et al., 1995, for a detailed description). The majority of the detected galaxies in the Survey were observed at multiple positions along the major axis, and the total CO fluxes were derived by modeling the underlying CO distribution.

The Survey is not complete in terms of any of the sample selection criteria, and the full sample cannot be used directly for a statistical analysis in a straightforward way. The analysis of the velocity integrated CO luminosity is further complicated by the fact that the CO line width varies from one galaxy to another (and even one position to another in a given galaxy) while the sensitivity achieved also varies depending on the observing condition even if identical integration times are used (typically 2-4 hours). Therefore, it is difficult to address the sensitivity and completeness limit of the Survey in terms of the observed CO flux. Instead, taking advantage of the well-known tight correlation between the FIR and CO luminosity (see Young et al., 1995), we utilize the IRAS 60 $`\mu `$m flux density of the individual sample galaxies to define a complete sample for a statistical analysis and the derivation of the CO LF. We further test the robustness of the technique and the selection bias by deriving the CO luminosity function using the optical $`B`$-band selection (see § 3.3). Previously Briggs & Rao (1993) successfully derived the HI mass function for the field galaxies using the optical selection, and the radio luminosity function has been successfully derived using the FIR selection function by Yun et al. (2001).

Taking advantage of the well defined and well studied IRAS Bright Galaxy Samples (BGSs; Soifer et al., 1989; Sanders et al., 1995), we adopt the sample selection criterion of 60 $`\mu `$m flux density limit greater than 5.24 Jy. A total of 200 galaxies in the Survey satisfy this criteria. Most of the selected galaxies are spirals of types Sa-Sc while 14 galaxies are mergers, 7 galaxies are close pairs, and 3 galaxies have no type determined. In most cases we adopt the 60 $`\mu `$m fluxes from the BGS surveys. For galaxies with angular diameters larger than 8′, we adopt the values from Rice et. al (1988) after multiplying by 1.18 in order to match the flux scaling (see Devereux & Young, 1990). This rescaling has little impact on our sample completeness limit because most of our sample galaxies have flux densities much larger than 5.24 Jy. Some of the 60 $`\mu `$m flux density measurements come from Young et al. (1989, 2002), and 10 galaxies satisfying our selection criterion are added to the sample. Five of these galaxies are not in the BGS survey area, and two are present in Rice et. al (1988) but not in the BGS surveys.

There are 12 galaxies in in our selected sample that are brighter than 5.24 Jy at 60 $`\mu `$m and were not detected in CO (see Table 1). The CO (1–0) flux measurements for five of these galaxies are found in the literature. The remaining 7 galaxies are treated as two limiting cases: (a) as detections with zero flux; and (b) as detections at the upper limit flux value (see Verter, 1987). A large number of non-detections would severely limit the determination of the the true CO LF, but the non-detections account for less than 4% in our sample and thus have little overall impact.

Some of the galaxy distances come from direct measurements such as using Cepheids. We adopt a Hubble constant of $`H_0=75`$ km s-1 Mpc-1 for the remaining galaxies. Since we are using the CO database obtained from largely a single instrument (2.5% of the sample taken from other surveys), the internal consistency of the data and the analysis should be quite good. Our sample also includes 26 galaxies in the Virgo cluster. We adopt a uniform distance of 16 Mpc for the Virgo galaxies. We have constructed a CO LF with and without Virgo galaxies in order to see if the presence of the cluster galaxies influenced the shape of the CO LF (see § 3).

To examine the uniformity and completeness of our sample, we have analyzed the differential number count statistics as a function of the $`60\mu m`$ flux density (see Fig. 1). If the $`60\mu m`$ sources are uniformly distributed in a Euclidean space and are not evolving, the resulting differential number count should be a power-law distribution with $`NS^{3/2}`$. The number of sample galaxies per bin shown in Figure 1 is consistent with such a power-law, but the slope is somewhat shallower. When all of the flux density bins are included, the best fit power-law index is $`\alpha =1.0`$, where $`logN=N_0+\alpha \times logS`$. For the $`60\mu m`$ flux density between 10 and 100 Jy where most of the sample galaxies fall, the differential number counts are more consistent with the uniform distribution within the statistical uncertainties.

There are several factors contributing to the flatter than expected power-law slope in the observed number counts in Figure 1. First of all, the spatial distribution of galaxies is in general not uniform, particularly given the relatively high flux density cutoff we adopted. The situation is exacerbated by the fact that our Galaxy is located within a galaxy aggregate called the Local Supercluster. The flattening of the galaxy counts among the high flux density bins as seen in Figure 1 is an immediate outcome of the local large scale structures. Soifer et al. (1989) noted a similar number count enhancement in their analysis of the first Bright Galaxy Sample (BGS1). The mean $`<V/V_m>`$ for our sample is about 0.35, indicating that on average our sample galaxies are about 10% closer than the uniform case. The $`<V/V_m>`$ values are particularly low for the log $`L_{60\mu m}`$ bins of 9.6-10.4 (see Fig. 2), and similar trend is also seen in the $`<V/V_m>`$ plot by Soifer et al. (1989). Comparison to Soifer et al. (1989) value $`<V/V_m>=0.47`$ is suggesting that major cause of this discrepancy is weak bias in our sample towards galaxies closer than average expected from the complete sample, while the large scale structure plays less important role. Although the FCRAO CO Survey selected target galaxies “at random” from the parent sample of IRAS and $`B`$-band selected galaxies, a slight bias favoring galaxies which were well suited to the FCRAO resolution (1-5 in size) is also present. Such a bias manifests non-trivially in these statistics, however. Given the limitation of our sample, such as selecting galaxies suitable for observation with FCRAO telescope and the local large scale structures, the steep power-law index of $`\alpha 1.2`$ shown in Figure 1 suggests that our assumption of random and uniform sampling from the parent samples still seems reasonable. Effects of these potential biases on the derivation of the LFs are evaluated by the derivation of the 60$`\mu m`$ LF using the same sample, as discussed bellow.

## 3 Luminosity Functions

A luminosity function represents the space density of the sample galaxies per luminosity bin ($`\mathrm{\Delta }L`$) centered on L. It represents the probability density of finding a galaxy with a specific luminosity, and it also contains information on the total luminosity in the sampled volume. We use the classical $`1/V_m`$ method (Schmidt, 1968) to derive the LF and associated uncertainty directly as:

$$\rho (L)=\underset{i=1}{\overset{N}{}}\frac{1}{V_m^i},\sigma _\rho (L)=(\underset{i=1}{\overset{N}{}}\frac{1}{(V_m^i)^2})^{1/2}$$
(1)
where $`\rho (L)`$ is number of objects per volume per luminosity bin centered on luminosity L, $`V_m`$ is the sample volume appropriate for each of the galaxies, and N is the number of galaxies in the bin.

The FCRAO Extragalactic CO Survey covers the same area as the IRAS BGS (Soifer et al., 1989; Sanders et al., 1995), and we selected a subset using the same IRAS 60 $`\mu `$m flux density cutoff. Therefore, the sampling correction for our sample can be reduced to a simple scaling relation if our sample is a fair subset of the BGS, as argued above (§ 2). In this case, the scaling factor is the ratio of the number of galaxies contained in BGS surveys and number of galaxies in our sample.

The sample volume obtained using this approach is:

$$V_m=\mathrm{\Omega }\frac{N_{FCRAO}}{N_{BGS}}\frac{D^3}{3}=1.16\times D^3$$
(2)
where D is maximum distance at which galaxy of given luminosity could be detected for a given flux limit, and $`\mathrm{\Omega }`$ is the total solid angle covered by both BGS1 (Soifer et al., 1989) and BGS2 (Sanders et al., 1995) in sterradian. $`N_{FCRAO}/N_{BGS}`$ is the ratio of the number of galaxies in our sample and the number of galaxies in both BGS samples. BGS1 and BGS2 together contain 601 galaxies and cover 83% of the sky.

### 3.1 60 $`\mu `$m Luminosity Function

The 60 $`\mu `$m luminosity represents the contribution from the IRAS 60 $`\mu `$m band to the FIR luminosity (Helou et al., 1988) and was calculated using the following relation from Yun et al. (2001) :

$$logL_{60\mu m}(L_{})=6.014+2logD+logS_{60\mu m}$$
(3)
where D is the distance in Mpc and $`S_{60\mu m}`$ is the IRAS 60 $`\mu `$m band flux density in Jy.

As a way to verify our CO sample correction and to evaluate the effects of various sample biases discussed in § 2, we first derive the 60 $`\mu `$m LF for our CO-selected galaxies in order to compare it to the previously published 60 $`\mu `$m LF by Yun et al. (2001). As shown on Figure 3, the agreement between the two 60 $`\mu `$m luminosity functions is excellent. Yun et al. sample was selected with 2.5 times smaller flux density cutoff than the BGS and is thus 9 times larger in size. While this is not a complete verification of our derivation method, this agreement lends support to our assumption that our sample is a fair subset of the complete 60$`\mu m`$ selected sample. The normalization factor could in principle be dependent on the 60$`\mu m`$ luminosity. The comparison of the 60$`\mu m`$ LFs suggests that this dependence is weak if any.

About 13% of our sample are Virgo cluster galaxies. Again, the 60 $`\mu `$m luminosity function was constructed including and excluding Virgo galaxies. The resulting luminosity functions are very similar – see Figure 4. The similarity would arise if the cluster galaxies are distributed in luminosity bins in the same manner as the field galaxies and if the Virgo galaxies contribute only a small fraction to each luminosity bin. Previous studies of Virgo galaxies have suggested that their molecular gas contents are affected little by their cluster environment (e.g. Kenney & Young, 1989). Since we are interested in the total cold gas density of the local volume regardless of the galaxy environment, we include all Virgo galaxies in our analysis here on.

### 3.2 CO Luminosity Function

Encouraged by the success of deriving the 60 $`\mu `$m LF above for the 200 galaxies in our sample, the CO luminosity function is derived using the 12CO(1–0) flux with the same sample and volume correction. The CO luminosity is calculated as $`L_{CO}=4\pi D^2S_{CO}`$, where D is the distance to the galaxy in Mpc and $`S_{CO}`$111We adopt the definition $`S_{CO}I_{CO}𝑑\mathrm{\Omega }`$, following the same notation as Young et al. (1995). is the velocity integrated CO flux in Jy km s-1. The computed CO luminosity for our sample galaxies ranges between $`10^4`$ and $`10^8`$ Jy km s-1 Mpc2. They are divided into 10 bins of unit magnitude (0.4 wide in a log scale). There are larger uncertainties at the low luminosity end because only a small number of sources occupy these bins.

The resulting CO luminosity function is shown in Figure 4. It is well described by a Schechter function which has an exponential cutoff at high luminosity end. The low luminosity end power-law slope is quite flat, similar to the 60 $`\mu `$m LF (see § 4). To quantify these properties, the characteristic parameters are derived for the Schechter function (Schechter, 1976) of the form

$$\varphi (L)d(L)=\varphi ^{}(\frac{L}{L^{}})^\alpha exp(\frac{L}{L^{}})d(\frac{L}{L^{}})$$
(4)
where $`\varphi ^{}`$ is the normalization factor, $`\alpha `$ is the low luminosity end power-law slope, and $`L^{}`$ is the characteristic luminosity at which the exponential cutoff develops. Since we are using logarithmic intervals to construct LF, we need to change this equation to an appropriate form:

$$\rho (L)=\rho ^{}(\frac{L}{L^{}})^{\alpha +1}exp(\frac{L}{L^{}})\mathrm{ln}10$$
(5)
in this case $`\rho (L)=\varphi (L)d(L)/d(\mathrm{log}(L))`$, and it is calculated using eq. ( 1). The derived best fit Schechter parameters are: $`\rho ^{}=(0.00072\pm 0.00035)`$ Mpc-3 mag-1, $`\alpha =(1.30\pm 0.16)`$, and $`L^{}=(1.0\pm 0.2)\times 10^7`$ Jy km s-1 Mpc2 ($`M^{}=9.4\pm 1.9\times 10^9M_{}`$ using the standard CO-to-H2 conversion factor – see § 4.2).

The single Schechter function fit using all 10 luminosity bins produces a rather poor fit with $`\chi ^2=50`$. This poor fit is also obvious in Figure 4 as the single Schechter function fit shown with a solid line is only marginally consistent with the data points and their formal uncertainties. It is possible that a Schechter function is intrinsically a poor functional form for the CO LF. On the other hand, a Schechter function successfully describes many other LFs including the HI LF (e.g. Schneider et. al, 1998). There is also a theoretical basis in that it is the functional form associated with the halo mass function in the Press-Schechter formulation of the structure formation (Press & Schechter, 1974).

One possible explanation for such a large $`\chi ^2`$ value is that there are two (or more) distinct populations of galaxies in the sample, similar to the situation with the FIR luminosity function as noted by Yun et al. (2001). In the FIR, the light-to-mass ratio for the starburst population is systematically enhanced by 1-2 orders of magnitudes over the field population. The CO luminosity is also thought to be elevated significantly among the intense nuclear starburst systems (Scoville, Yun, & Bryant, 1997; Downes & Solomon, 1998), and this may materialize as two distinct populations in the CO LF as well. There are 14 mergers/starburst systems in our sample, and eight of these fall within the last two bins and account for about 40% of the galaxies in these bins. When a best fit Schechter function is derived using only the first 8 bins dominated by normal field galaxies, the $`\chi ^2`$ value drops to 13, which is a significant improvement. The Schechter function fit for the bins dominated by normal galaxies gives: $`\rho ^{}=(0.003\pm 0.001)`$ Mpc-3 mag-1, $`\alpha =(0.90\pm 0.15)`$, and $`L^{}=(3.3\pm 0.6)\times 10^6`$ Jy km s-1 Mpc2 ($`M^{}=3.1\pm 0.6\times 10^9M_{}`$).

The difference between all Schechter parameters for the two limiting situations of treating the non-detected galaxies as upper limits or zero detections is smaller than 2%. Such a small difference is expected since we have less than 4% upper limits in our sample. Similarly, including or excluding Virgo galaxies has only a minor impact, as discussed above (see § 3.1).

### 3.3 CO LF from the $`B`$-band Selected Galaxies

As an independent check of our CO LF derived using the FIR selection, we also derive a CO luminosity function using an optical $`B`$-band selected sample. The method of deriving the CO LF is similar to the method used for the infrared selected galaxies: we use the $`1/V_m`$ method and rescaling of the number of galaxies using a complete sample, which is the Revised Shapley-Ames (RSA) catalog by Sandage & Tammann (1981) in this case. A declination cut off of $`\delta >25^{}`$ and a Galactic longitude limit of $`|b|>10^{}`$ are adopted as part of the sample definition. This makes the total surveyed area equal to 60% of the sky. Our brightness selection criterion is $`B_T12.0^{mag}`$.

Accounting for different Hubble types correctly is a difficult issue. Since CO emission is generally very weak or undetected among early type galaxies, only the spiral and irregular galaxies are included in the analysis. The examination of the differential number counts for our $`B`$-band selected sample yields a power-law index $`\alpha =1.35`$, suggesting a reasonably good uniformity and completeness. The mean $`<V/V_m>`$ ratio for the sample is 0.33 (see Fig. 2), similar to the FIR selected sample. Although we adopt the RSA catalog as a complete sample, a few galaxies satisfying our selection criteria are missing in the RSA (see Condon, 1987).

The $`B`$-band magnitudes are taken from the RSA while the CO fluxes are taken from the Survey. There are 133 galaxies in the CO Survey that satisfy the selection criteria while the parent RSA catalog galaxies satisfying the selection criteria are 252 in total. In this case, CO fluxes for 13 galaxies (10%) are upper limits. The CO luminosity function derived for this $`B`$-band selected sample is shown on the Figure 5. This LF is very similar to CO LF from the FIR selected galaxies but with a lower characteristic luminosity $`L^{}`$ and a larger scatter on the low luminosity end. The lower characteristic luminosity $`L^{}`$ is a direct consequence of the absence of mergers in this sample. The best fit Schechter parameters for the $`B`$-band derived CO luminosity function are: $`\rho ^{}=(0.0021\pm 0.0009)`$ Mpc-3 mag-1, $`\alpha =(1.0\pm 0.2)`$, and $`L^{}=(4.8\pm 1.1)\times 10^6`$ Jy km s-1 Mpc2 ($`M^{}=4.5\pm 1.0\times 10^9M_{}`$).

## 4 Discussion

### 4.1 Uncertainties and biases in the CO LF

One major source of uncertainty in the derived CO luminosity function is that it is derived indirectly using the FIR and optical $`B`$-band selection functions. A direct derivation from the observed CO properties is possible in principle, but the complication associated with the a priori unknown line widths adds a significant uncertainty. Instead we take advantage of the known tight correlation between FIR and CO luminosity in deriving the CO LF for our sample. The 1.4 GHz radio luminosity function derived by Yun et al. (2001) using the radio-FIR correlation and IRAS 60 $`\mu `$m flux density agrees very well with the radio LFs derived directly (Condon, Anderson, & Helou, 1991; Condon et al., 2002), giving some assurance to this technique. We have examined whether any systematic trends (thus a bias) exist in the $`S_{60\mu m}/S_{CO}`$ ratio as a function of the 60 $`\mu `$m flux density and luminosity. As shown in Figure 6, little trend is seen in this ratio, and the known linear correlation seems to hold well within the uncertainties for the entire range of flux density and luminosity.

Deriving a luminosity function for one wavelength using a selection function at another wavelength may be robust enough to work even if the two quantities are correlated in a non-linear way. There is a known correlation between FIR, CO, and optical $`B`$-band luminosity for late type galaxies, but the correlations involving the $`B`$-band luminosity are shown to be non-linear (see Young et al., 1989; Perea et al., 1997). Yet, the CO luminosity function derived using the FIR selection function (Figure 4) is in excellent agreement with the CO LF derived using the $`B`$-band selection function (Figure 5). The selection bias is present as some of the most CO luminous galaxies (merger starbursts) are missed in the $`B`$-band selected sample. Nevertheless the agreement is striking given the non-linear dependence between CO and $`B`$-band luminosity.

Another source of a significant uncertainty is the uniformity of sampling and correctly accounting for the selection function for the sample of galaxies which is not complete by the adopted selection criteria. Both the differential source count and the $`<V/V_m>`$ analysis suggest a non-uniform distribution of the sample galaxies (see § 2), and a slight bias towards galaxies brighter at 60$`\mu m`$ may be responsible for this effect. The large scale structure also have influence in lowering $`<V/V_m>`$ value for our sample. The non-trivial nature of these effects are demonstrated by the fact that including and excluding Virgo cluster galaxies from the sample make little difference to the differential source count and the $`<V/V_m>`$ values.

### 4.2 Local Cold Gas Mass Density

CO emission is a commonly used tracer of molecular hydrogen because CO is one of the most abundant molecules in cold ISM and because its excitation in astrophysical conditions is determined by collisions with hydrogen molecules. Using the CO LF derived above, we can for the first time estimate the total mass density of molecular gas in the local volume. Adopting $`N(H_2)/I(CO)=3\times 10^{20}`$ cm-2 \[K km s-1\]-1 (see review by Young & Scoville, 1991),

$$M(H_2)=1.18\times 10^4S_{CO}D^2M_{}$$
(6)
where $`S_{CO}`$ is total CO flux of a galaxy in Jy km s-1 and $`D`$ is luminosity distance in Mpc.

The molecular mass density in the local volume contributed by each luminosity bin is shown in Figure 7. The dominant contribution to the mass density comes from galaxies around $`L^{}`$ as expected. The summation over the 10 bins gives $`\rho _{H_2}=M_{H_2}/V_m=(2.4\pm 0.7)\times 10^7M_{}`$ Mpc-3. Integration of the LF using the Schechter parameters obtained in previous section gives a little bit smaller values: $`\rho _{H_2}=(2.2\pm 1.1)\times 10^7M_{}`$ Mpc-3 for the fit trough all 10 bins, and $`\rho _{H_2}=(2.2\pm 0.9)\times 10^7M_{}`$ Mpc-3 for fit trough first 8 bins. The second fit is more realistic, since it fits much better bins that dominate contribution to the total mass, i.e. bins around $`L^{}`$. The uncertainty stated is $`\pm 1\sigma `$. The systematic uncertainty, which includes uncertainties in the flux measurements, distance determinations, and CO-to-H2 conversion, is probably larger. Unless otherwise is stated, we adopt $`\rho _{H_2}=(2.3\pm 0.9)\times 10^7M_{}`$ Mpc-3 as an average value between values obtained from the fit and the direct summation. Since the dependence of the gas mass density on the Hubble constant ($`hH_0/100`$ \[km s-1 Mpc-1\]-1) is linear, this result can be written as $`\rho _{H_2}=(3.1\pm 1.2)\times 10^7hM_{}`$ Mpc-3.

The molecular hydrogen mass density in the local volume obtained from the $`B`$-band selected sample is $`\rho _{H_2}=(3.1\pm 0.9)\times 10^7hM_{}`$ Mpc-3 using a direct summation of the contribution of each galaxy and $`(3.1\pm 1.5)\times 10^7hM_{}`$ Mpc-3 using the best fit Schechter parameters. Values for the local molecular mass density obtained from the FIR selected sample and the $`B`$-band selected sample are in good agreement. Since we have better statistics for the larger FIR selected sample, we adopt the local molecular gas mass density derived from the FIR selected sample here on. Using observed H2/HI mass ratio for different morphological types of galaxies and fraction of each morphological type Fukugita et al. (1998) estimated a similar value for $`\rho _{H2}`$.

In comparison, Zwaan et al. (1997) estimated the local atomic gas mass density of $`\rho _{HI}=(5.8\pm 1.2)\times 10^7hM_{}`$ Mpc-3 from their Arecibo HI strip survey while Rao & Briggs (1993) derived $`\rho _{HI}=(4.8\pm 1.1)\times 10^7hM_{}`$ Mpc-3 from a sample of optically selected galaxies. Therefore, the molecular gas mass density in the local volume is about 50-65% of the atomic mass density, and molecular gas represents a significant component of the total mass density of the neutral gas. The HI masses for 176 galaxies in our sample are known (Young et al., 1995), and we have computed the HI mass density in the local volume using the same procedure as for the $`H_2`$ mass density. A summation over the HI mass weighted by $`V_m`$ gives a value consistent with the value obtain by Zwaan et al. (1997).

We examined our sample for any trends in $`M_{H_2}/M_{HI}`$ ratio vs. $`M_{HI}`$ and $`M_{H_2}`$. For the range of HI masses, $`8.8<logM_{HI}(M_{})<11`$, the ratio of molecular to neutral atomic hydrogen masses is on average around unity with considerable scatter (see Fig. 8a). Masked in the large scatter is a possible trend in $`M_{H_2}/M_{HI}`$ as a function of Hubble type (see Young et al., 1989). A clearer trend of increased $`M_{H_2}/M_{HI}`$ ratio with increasing H2 mass is seen in Figure 8b. For $`logM_{H_2}(M_{})<8.5`$, the average of this ratio is below 0.5 while the ratio jumps to around 2 near $`M_{H_2}10^9M_{}`$. Above H2 mass of $`10^9M_{}`$, on average galaxies have more molecular than atomic gas.

Using the molecular gas mass density derived here and the value of HI mass density from Zwaan et al. (1997) we derive a total cold gas mass density at present epoch as a fraction of the critical density $`\mathrm{\Omega }_{HI+H_2}=(3.2\pm 0.8)\times 10^4h^1`$. Supposing that He contributes 25% of the total gas mass density, we add 33% to the derived value of $`\mathrm{\Omega }_H`$, which gives $`\mathrm{\Omega }_{HI+H_2+He}\mathrm{\Omega }_{gas}=(4.3\pm 1.1)\times 10^4h^1`$. These values for total neutral gas mass density are 15% smaller if we use Rao & Briggs (1993) value for the HI mass density.

The present baryon density estimated from the D/H ratio and Big Bang nucleosynthesis is $`\mathrm{\Omega }_b=0.02\pm 0.002h^2`$ (Burles et. al, 2001) while the estimate from the cosmic microwave background anisotropy is also around $`\mathrm{\Omega }_b0.02h^2`$ (de Bernardis et. al, 2002). The stellar mass contents account for $`\mathrm{\Omega }_{}0.0025h^1`$ (Cole et al. 2001; but also see Benson, Frenk, & Sharples 2002). We conclude that the cold gas content of late type galaxies is around 20% of the stellar mass content and about 2% of the total baryonic content in the current epoch.

### 4.3 CO-to-H2 conversion

The lack of electric dipole moment makes direct observations of molecular hydrogen difficult in general, and studying the spatial extent and molecular gas mass requires another tracer. Highly abundant, chemically robust, and easily excited by collision with H2 molecules, CO is the most commonly used tracer of molecular gas. The CO (1–0) transition is optically thick under most astrophysically interesting conditions, making it relatively insensitive to metallicity and abundance effects. The two key excitation parameters of density and temperature have a nearly canceling effect, making CO a fairly reliable tracer of H2 in a broad range of physical conditions (see reviews by Maloney & Black, 1988; Young & Scoville, 1991).

The derived ratios of $`N(H_2)/I(CO)`$ range between $`(15)\times 10^{20}`$ cm-2 \[K km s-1\]-1 (Bloemen et al., 1986; Dickman et al., 1986; Scoville & Sanders, 1987). We adopt a constant CO-to-H2 conversion factor of $`N(H2)/I(CO)=3\times 10^{20}`$ cm-2 \[K km s-1\]-1 (see discussions by Young & Scoville, 1991). Devereux & Young (1990) show that the H2 mass estimates of galaxies from their CO luminosity are accurate to $`\pm 30\%`$. Young & Scoville (1991) show that the CO-to-H2 conversion for galaxies of diverse morphology and metallicity are similar in absolute value to the conversion in the Milky Way. While the $`H_2`$ mass estimate for an individual galaxy may be uncertain to about 30%, the $`H_2`$ mass estimate for an ensemble of galaxies should be more reliable. In most galaxies, giant molecular clouds (GMCs) and cloud complexes dominate the total molecular gas mass, and adopting a conversion factor consistent with the values derived from the Galactic GMCs should yield the most robust, mass-weighted estimates of gas masses.

Among the low metallicity galaxies such as Small Magellanic Cloud, CO abundance may become low enough to affect the standard assumption of self-shielding and thermalization. The derived conversion factors are generally larger (e.g. Wilson, 1995), and low metallicity dwarf galaxies are often undetected entirely in CO. Therefore the H2 mass and mass density derived from the CO luminosity are strictly lower limits since molecular gas from these galaxies is missing from our analysis. On the other hand, the contribution by the low luminosity bins to our derived $`H_2`$ mass density is small (see Fig. 7) as the majority of the contribution to the total mass density comes from bins around $`L^{}`$. Therefore, we can largely neglect the effect of low metallicity, low luminosity galaxies in the derivation of the total molecular gas mass density for the local volume.

As discussed already in § 3.2, the CO emission may be elevated for the FIR luminous galaxies with nuclear gas concentrations, and the CO-to-H2 conversion factor may be smaller than the canonical value. In the extreme environment of the nuclear starburst regions, CO emission arises from a multi-phase medium with sub-thermal excitation in the diffuse phase, and the standard conversion factor may over-estimate the molecular gas mass by a factor as large as 3 to 5 (Scoville, Yun, & Bryant, 1997; Downes & Solomon, 1998). Such merger/starburst galaxies contribute significantly only for the top two luminosity bins in our CO LF, and these two bins contribute less than 5% to the the total molecular gas mass density. Therefore, the high CO luminosity conversion factor does not have a significant impact on the derived total molecular mass density.

Similarly, increased CO emission in the central 1 kpc of our Galaxy and other galaxies has been suggested by recent observations (Paglione et. al, 2001, S. Hüttemeister, private communication). High angular resolution observations of gas-rich spiral galaxies frequently reveal a distinct component in the central 1 kpc (see Sakamoto et al., 1999), and some fraction of the total CO luminosity may arise from such a component. On the other hand, among the Virgo spirals studied by Kenney & Young (1988), a simple exponential distribution without a significant central component offers a good fit for the 12 out of 14 galaxies whose CO emission is spatially well resolved (i.e., CO detected at $`4`$ positions along the major axis). We made no attempt to account for enhanced CO emission in the central kpc as the required information is generally not available. If the molecular gas properties of Virgo spirals are typical of the late type field galaxies, then the possible contribution by the enhanced nuclear CO emission in some galaxies may not be substantial.

In summary, using the canonical CO-to-H2 conversion relation is problematic in some cases, such as low metallicity systems or luminous nuclear starburst systems, and CO is a poor tracer of molecular gas among low metallicity dwarfs. For these reasons the local molecular gas mass density we infer is really a lower limit. However, since the majority of the mass contribution comes from $`L^{}`$ galaxies whose molecular gas mass is dominated by GMCs like our Galaxy, the use of the standard conversion factor still offers a fairly reliable estimate of the total molecular gas density for the local volume.

## 5 Summary

Utilizing the largest available CO survey, the FCRAO Extragalactic CO Survey (Young et al., 1995), the CO luminosity function for the local volume is derived using (1) a FIR selected sample of 200 galaxies that satisfy $`S_{60\mu m}>5.24Jy`$ and (2) optical $`B`$-band selected sample of 133 galaxies. Although neither of the samples is complete in terms of the sample selection, a sampling function is constructed using a well defined parent sample, and a non-parametric CO luminosity function is derived from each sample. By examining the properties of the CO luminosity functions, we conclude:

1. The CO luminosity functions derived from the FIR and $`B`$-band selected samples are reasonably well described by a Schechter function. The characteristic luminosity $`L^{}`$ is around $`L_{CO}10^7`$ Jy km s-1 Mpc2. The low luminosity end of the CO luminosity function is ($`\alpha =1.3\mathrm{to}0.9`$). Similar values are obtained for the two CO LFs derived using the two differently selected samples.
2. The molecular gas mass density of the local volume is $`\rho _{H_2}=(3.1\pm 1.1)\times 10^7hM_{}`$ Mpc-3 which is about 50-65% of the HI gas mass density. This value is not strongly affected by variations in the CO-to-H2 conversion factor for low and high luminosity galaxies since it is dominated by the $`L^{}`$ galaxies with $`L_{CO}10^7`$ Jy km s-1 Mpc2 ($`M_{H_2}^{}5\times 10^9M_{}`$). Its dependence on the global CO-to-H2 conversion factor is linear and should be secure to within a factor of 30% or better.
3. Combined with the HI gas mass density, we estimate the total cold gas mass density at the present epoch as a fraction of the critical density $`\mathrm{\Omega }_{HI+H_2}=(3.2\pm 0.8)\times 10^4h^1`$. When the He contribution is included, the cold gas mass density increase to $`\mathrm{\Omega }_{HI+H_2+He}\mathrm{\Omega }_{gas}=(4.3\pm 1.1)\times 10^4h^1`$. Therefore, the cold gas content of late type galaxies corresponds to about 20% of the stellar mass content and about 2% of the total baryonic content in the universe.

The authors acknowledge insightful discussions with N. Katz, M. Heyer, S. Hüttemeister, and others. Comments by our referee, Martin Zwaan, were quite helpful for improving the manuscript. This research has made use of the NASA/IPAC Extragalactic Database (NED) which is operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with the National Aeronautics and Space Administration. D. Keres gratefully acknowledge the partial research support provided by the Mary Dailey Irvine graduate research fellowship. The Five College Radio Astronomy Observatory is operated with the permission of the Metropolitan District Commission, Commonwealth of Massachusetts, and with the support of the National Science Foundation under grant AST 97-25951.