{'fraction_non_alphanumeric': 0.040450257787503524, 'fraction_numerical': 0.02292733873775302, 'mean_word_length': 4.435393764642278, 'pattern_counts': {'":': 0, '<': 6, '<?xml version=': 0, '>': 14, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Near Infrared Imaging Survey of Bok Globules: Density Structure

On the basis of near-infrared imaging observations, we derived visual extinction (Av) distribution toward ten Bok globules through measurements of both the color excess (E_{H-K}) and the stellar density at J, H, and Ks (star count). Radial column density profiles for each globule were analyzed with the Bonnor-Ebert sphere model. Using the data of our ten globules and four globules in the literature, we investigated the stability of globules on the basis of xi_max, which characterizes the Bonnor-Ebert sphere as well as the stability of the equilibrium state against the gravitational collapse. We found that more than half of starless globules are located near the critical state (xi_max = 6.5 +/- 2). Thus, we suggest that a nearly critical Bonnor-Ebert sphere characterizes the typical density structure of starless globules. Remaining starless globules show clearly unstable states (xi_max>10). Since unstable equilibrium states are not long maintained, we expect that these globules are on the way to gravitational collapse or that they are stabilized by non-thermal support. It was also found that all the star-forming globules show unstable solutions of xi_max>10, which is consistent with the fact that they have started gravitational collapse. We investigated the evolution of a collapsing gas sphere whose initial condition is a nearly critical Bonnor-Ebert sphere. We found that the column density profiles of the collapsing sphere mimic those of the static Bonnor-Ebert spheres in unstable equilibrium. The collapsing gas sphere resembles marginally unstable Bonnor-Ebert spheres for a long time. We found that the frequency distribution of xi_max for the observed starless globules is consistent with that from model calculations of the collapsing sphere.


Introduction
Probing the physical evolution of dense cores toward the onset of star formation is of great importance in understanding the nature of the star formation process. The star formation process as well as the characteristics of forming stars are probably mostly determined by the properties of the nursing cores, such as size, mass, temperature, turbulence, and density structure, right before the gravitational collapse. In order to clarify the initial condition of star formation, it is important to provide a number of dense core samples with well defined physical properties.
Bok globules (i.e., isolated molecular cloud cores) are ideal sites for studying isolated star formation processes because of their simple shape and isolated geometry from neighboring clouds. These characteristics of globules enable us to model their internal structure (e.g., density). For these reasons, detailed observations of some well-known globules have been done. Alves et al. (2001) reported a high-resolution dust-extinction (A V ) study of the starless globule Barnard 68, based on near-infrared observations. They presented a radial column density profile of the globule which is well fitted by the Bonnor-Ebert sphere model (Bonnor 1956;Ebert 1955). Because the density structure is well determined, molecular abundances in Barnard 68 have been studied by comparing dust extinction with molecular distribution Hotzel et al. 2002aHotzel et al. , 2002bDi Francesco et al. 2002;Lai et al. 2003).
An unbiased survey of both starless and star-forming globules is important to provide essential information on the evolution of density, velocity, and chemical structure of globules as well as their relative lifetime. Statistical studies of globules (dense cores) have been conducted based on the observations of thermal dust emission (e.g., Visser et al. 2001Visser et al. , 2002, dust extinction in the mid-infrared (Bacmann et al. 2000), molecular line emission (e.g., Launhardt et al. 1998;, and a combination of these (e.g., Tafalla et al. 2002; dust and molecular line emissions). However, no statistical studies of the density structure of globules based on the near-infrared extinction have been carried out. Recently, the dust temperature gradient in cold starless cores (e.g., Evans et al. 2001;Ward-Thompson, Andre, & Kirk 2002) has been pointed out, and this may systematically affect the column density distribution derived by assuming a constant dust temperature. The measurements of extinction in the mid-infrared surface photometry (e.g., Bacmann et al. 2000) can also be affected by dust temperature distributions, because both absorption and emission from the line-of-sight medium should be taken into account at these wavelengths. As a column density tracer, near-infrared extinction measurement is one of the most straightforward methods with no ambiguity caused by dust temperature distribution. Near-infrared extinction measurements can trace dust column densities at low-density outer regions of the globules, where (sub)millimeter dust emission is hard to detect.
The aim of our near-infrared (J, H, and K s ) imaging survey of Bok globules is to reveal the evolution of radial density structures of globules from the starless to protostellar phases as well as to derive their fundamental physical parameters (e.g., size, temperature, mass, and external pressure) on the basis of the measurements of dust extinction. We investigated the radial A V distribution for each globule using the Bonnor-Ebert sphere model. The solutions of the Bonnor-Ebert sphere can be characterized by the dimensionless radial parameter ξ max , which describes the stability of the gas sphere against gravitational collapse, and the critical state is achieved at ξ max =6.5. By -4combining our results for ten Bok globules with previous reports on four sources, we investigate the difference in physical properties, in particular the ξ max value, between the globules in the starless phase and those in the star-forming (protostellar or gravitational collase) phase, which is likely to be related to the difference in evolutionary stages of the globules. In order to investigate the velocity structure of globules, we observed the same globules using the 45 m radio telescope of Nobeyama Radio Observatory with the C 18 O (J = 1 → 0) and/or N 2 H + (J = 1 → 0) molecular lines.
The procedures of our observations and data reduction are described in § 2. The distribution of dust extinction (A V ) for each globule is obtained in order to derive the radial A V profile and to fit it by the Bonnor-Ebert sphere model ( § 3 and 4). We investigate the difference in the physical properties between the starless and star-forming globules and discuss implications for the evolution of globules in § 5. Our conclusions are summarized in § 6.

Source Selection
We observed ten Bok globules, CB 87,CB 110,CB 131,CB 134,CB 161,CB 184,CB 188,Lynds 495, and Lynds 498. These were selected from a unified catalog of dark clouds (Dutra & Bica 2002; compilation of 21 dark cloud catalogs) with the following criteria: (1) globules whose angular sizes are less than ∼5 ′ , (2) those having simple shapes (circular or elliptical), (3) those well isolated from neighboring clouds, (4) those located in the solar neighborhood (distance of less than ∼500 pc), and (5) those observable from both hemispheres. Our selection criteria of (1) to (4) are similar to those used in the CB catalog (Clemens & Barvainis 1988), which lists 248 small, isolated, optically selected dark clouds. Typical angular size and ellipticity of dark clouds in the CB catalog are ∼ 3 ′ and ∼ 2, respectively (see, Figures 1 and 2 of Clemens & Barvainis 1988). In the CB catalog, the fraction of large (e.g., ≫ 5 ′ ) and/or highly elliptical (e.g., ellipticity ≫ 2) clouds with isolated geometry is small. The properties of the observed globules are summarized in Table 1. Column 1 lists the name of the globule, Columns 2 and 3 list the equatorial coordinates (J2000) of the center of the target field. Column 4 gives the distance to the globule from the catalog (Dutra & Bica 2002). The distance to the nearby and compact globules is difficult to determine. The lack of foreground stars toward globules makes spectro-photometric distance measurements impossible, and the distance estimation from LSR velocity (kinematic distance) is not reliable for nearby clouds. Thus, the listed distances are only values assumed by associating observed globules with neighboring larger molecular clouds whose distance is known. Column 5 represents the association of IRAS point sources within the optical boundary of each globule. We classified globules into starless and star-forming clouds based on the association of IRAS point sources * . A -5 -Class I protostar (IRAS 19179+1129) is associated with CB 188 (Launhardt 1996;Yun & Clemens 1992), whereas no evidence of star formation is shown for the other nine globules. For CB 184, a Class II source (IRAS 19116+1623) is located near the globule (Launhardt 1996), but its physical association with the globule is not known. Since star-forming activity is not taking place around the center of CB 184, we will treat CB 184 as a starless globule in this paper. The Bonnor-Ebert model fitting studies of globules based on near-infrared extinction have been reported in four sources in the literature: Barnard 68 , Barnard 335 (Harvey et al. 2001), Coalsack Globule II (Lada et al. 2004;Racca et al. 2002), and Lynds 694-2 (Harvey et al. 2003). A Class 0 protostar (IRAS 19345+0727) is associated with Barnard 335. The other globules have no IRAS sources, while Lynds 694-2 shows strong evidence of gas inward motion (Lee, Myers, & Tafalla 1999Lee, Myers, & Plume 2004). Thus, we will classify Lynds 694-2 as a star-forming globule in this paper. There are two recent reports on the density structure of Coalsack Globule II based on similar analyses. We mainly refer to the observations of Lada et al. (2004) because their near-infrared observations are more sensitive than those of Racca et al. (2002).

Near-Infrared Data
We carried out near-infrared (J, H, and K s ) imaging observations of Bok globules using the infrared camera SIRIUS (Simultaneous three-color InfraRed Imager for Unbiased Surveys; Nagayama et al. 2003;Nagashima et al. 1999) on the 1.4 m telescope, IRSF (InfraRed Survey Facility), at the South African Astronomical Observatory (SAAO). Observations were made over two periods from 2002 July 30 to August 12 and from 2003 June 24 to July 14. SIRIUS has three 1024×1024 pixel HgCdTe infrared detectors (HAWAII array), which enables us to take J, H, and K s images simultaneously. IRSF/SIRIUS can image a large FOV of ∼7. ′ 7×7. ′ 7 with a scale of 0. ′′ 45/pixel, and the limiting magnitudes (S/N=10) for a typical 10-minute exposure time reach 18.9, 18.3, and 17.3 mag at J, H, and K s , respectively. Thus, this survey is much deeper than previous extensive surveys (e.g., 2MASS, DENIS). For each globule, we performed a 15-minute integration with a set of dithered 20 s or 30 s exposures except CB 110 and FeSt 1-457, whose integration times are 45 and 60 minutes, respectively. The dithered spatial interval was 30 ′′ or 40 ′′ . Since the apparent diameter of each target globule is less than 5 ′ , we covered its entire extent at one time. The typical seeing during the observations was 0. ′′ 9−1. ′′ 4 (FWHM), corresponding to 2−3 pixels. To correct detector flatness and pixel-to-pixel variation, we obtained twilight flat and dark frames at the beginning and end of the nights. We observed some sets of near-infrared standard stars listed by Persson et al. (1998) during the observations for photometric calibration.
We processed the observed data in a standard way of near-infrared image reduction (i.e., dark subtraction, flat-field correction, and median sky subtraction) with the IRAF † (Image Reduction and Analysis Facility) and the IDL (Interactive Data Language, by Research Systems, Inc.) routines. After subtraction of averaged dark frame, each object frame was divided by normalized flat on the 24 µm image, it is not known whether the source is a protostar associated with CB 131.
-6frames. For correcting thermal emission pattern and the fringe pattern from the atomospheric OH emission, we constructed a median sky frame to be subtracted from the object frame. Since all the observed targets are located in the very crowded field of the Milky Way, we did not adopt the self-sky, with which the residual signal of (bright) stars cannot be excluded on the resulting self-sky frame. Therefore, we constructed the sky frame for each band by taking the median of dithered frames toward the different several (four to eight) objects observed in the same night, which resulted in a sky frame with no stellar residuals. This "average" median sky does not reflect the actual sky for each object because the median was taken across several object frames. We compared the "average" median sky frame with the self-sky frame for each object in order to correct their different sky levels, and the final (level-matched) median sky frames were made for each object. Though the thermal emission pattern is adequately corrected in the final median sky frames, the OH fringe pattern is time-averaged. We confirmed that the effect of the time-variation of the OH fringe pattern is within acceptable levels by checking the sky noise level of the final combined object frames. After the sky subtraction and frame-registration, we combined all the dithered images for each object with a 3σ clipping. We show three-color (JHK s ) composite images of globules in Figure 1.
For the source detection and photometry, we applied the IRAF/DAOPHOT package (Stetson 1987) to the images. We measured the instrumental magnitudes of the detected stars having a peak intensity greater than 4.0 σ noise level above the sky background by fitting a point-spread function (PSF) and subtracting it from the image. We iteratively performed stellar detections and psf-fitting photometries to the source-subtracted image in order to collect overlooked sources in the previous procedure. We manually checked the results of the photometries and removed some false detections by visual inspection. We summarize the results of stellar detection in Table 2. For calculating the plate solution, we matched the pixel coordinates of a number of detected stars with the celestial coordinates of their counterparts in the 2MASS point source catalog (Cutri et al. 2003) and obtained the transformation coefficient using the IRAF/IMCOORDS package, which resulted in an acceptable positional error of less than 0. ′′ 5 (rms). All the stars detected in our observations were cataloged separately with respect to the filters (JHK s ). These stars of single-band detections were used to generate an extinction map at each band with a star count method. We applied the following equations for the color conversion from the IRSF to the standard CTIO/CIT photometric system: These were derived from measurements of red standard stars (see, Nakajima et al. 2004 for details).

Molecular Line Data
We have carried out molecular line observations toward ten globules using the 45 m radio telescope of Nobeyama Radio Observatory over four observing periods (2003 March 10-19, 2003 April 8-10, 2004 January 9-12, and 2004 February 3-10). The observed objects are the same as those -7observed by IRSF/SIRIUS. The receiver employed was the 25-element focal-plane array receiver BEARS consisting of double-sideband SIS mixers (Sunada et al. 2000). All the target globules in Table 1 are mapped in position-switching mode in C 18 O (J = 1 → 0) at 109.782173 GHz (JPL catalog; see, http://spec.jpl.nasa.gov/) and/or N 2 H + (J = 1 → 0) at 93.1737767 GHz (Caselli, Myers & Thaddeus 1995) on a 13. ′′ 7 grid spacing. The half-power beamwidth for each element beam was estimated to be ∼17. ′′ 8 at 93 GHz and ∼15. ′′ 1 at 109 GHz. The map coverage (∼4 ′ ×4 ′ ) is sufficiently large compared with the extinction feature of globules. As a back end, we used a digital autocorrelator with a 37.8 kHz frequency resolution (∼0.12 kms −1 at 93 GHz, and ∼0.10 kms −1 at 109 GHz). The intensity was calibrated using the standard chopper-wheel method, and is expressed in terms of the corrected antenna temperature T * A . To correct for the sideband ratio, we calibrated T * A using the intensity obtained with the single-sideband receiver S100. The scaling factor from T * A of BEARS to that of S100 for N 2 H + (J = 1 → 0) was measured by ourselves, and that for C 18 O (J = 1 → 0) was provided by the observatory. A typical value of the scaling factor for each BEARS element was ∼ 1.5. We estimated the main beam efficiency of S100 at 93 and 109 GHz to be 0.515 and 0.442, respectively, by interpolating values at 86, 100, and 115 GHz provided by the observatory. The pointing accuracy was typically better than 7 ′′ by observing SiO maser sources at 43 GHz every 1-1.5 hours during the observations.

Derivation of A V and Molecular Line Parameters
A V measurement is one of the most straightforward ways to probe the density structure of dark clouds. There are two standard ways to do extinction mapping; one is the star count method (e.g., Dickman 1978;Cambrésy 1999;Dobashi et al. 2005) based on the measurements of the stellar density, and the other is the near-infrared color excess (NICE) method (Lada et al. 1994) based on measurements of the stellar color excess (e.g., E H−K ). By considering the stellar density in H and K s at the densest region of globules, we decided which method is better for each target globule. We analyzed CB 131, CB 134 and CB 188 with the star count method and analyzed the other globules with the NICE method.

Star Count
On the basis of the stars detected in the J, H, and K s bands, we derived extinction maps (A J , A H , and A Ks ) for CB 131, CB 134, and CB 188 using the star count method (see, e.g., Kandori et al. 2003). A logarithmic cumulative stellar density log N measured at the magnitude m λ,0 can be converted to extinction A λ as where ∆α and ∆δ are the offset positions from the equatorial coordinate of image center, and f −1 is the inverse function of f = log N ref (m λ ), which is the logarithmic cumulative luminosity function (i.e., Wolf diagram; Wolf 1923) constructed in the extinction-free reference field (see Fig.  2). Though f is often assumed to be linear as f = a + bm λ , we fitted it with a 4th-order polynomial because the function shows a slightly curved shape.
-8 -To measure stellar density distribution, log N (∆α, ∆δ, m λ,0 ), we arranged circular cells (typically ∼30 ′′ in diameter) separated by 9 ′′ intervals on each observed field (∼7 ′ ×7 ′ area), and counted the number of stars falling in the cells up to the threshold magnitude m λ,0 (see Table 2 for the definition) which is well above the limiting magnitudes. The diameter of circular cells depends on the stellar density in each cloud field. We smoothed the stellar density map with a Gaussian filter (FWHM=13. ′′ 5) to reduce noise. In order to derive A J , A H , and A Ks , we numerically calculated f −1 in equation (4) after substituting the actually measured stellar density log N (∆α, ∆δ, m λ,0 ) into f −1 . In the resulting extinction maps, A J tends to be saturated at the densest region of the cloud, but sufficiently reveals the diffuse extinction feature around the cloud, and A Ks has the opposite tendency. In order to relieve saturation, we made a composite A V map by combining the A J , A H , and A Ks data with the following equation: where α and β are the weighting parameters defined as denotes the maximum value of A H in each analyzing field. The R λV denotes the standard conversion factor from A λ to A V in Rieke & Lebofsky (1985), which is 3.55 and 5.71 for J and H, respectively. For K s , we assumed A V = 9.44A Ks on the basis of equation (3) and the conversion factors in Rieke & Lebofsky (1985). These conversion factors depend on the assumption of dust optical properties, in particular R V . We discuss systematic uncertainties in the derivation of A V in the Appendix. The resulting A V maps for CB 131, CB 134, and CB 188 are shown in Figure 2. Typical uncertainty of A V in the reference field (A V ∼ 0 mag) is estimated to be ∼0.9 mag.

Near-Infrared Color Excess
On the basis of the stars detected both in the H, and K s bands, we derived extinction (A V ) maps for CB 87, CB 110, CB 161, CB 184, FeSt 1-457, Lynds 495, and Lynds 498 using the NICE method. If we assume that the population of stars toward the observed field is invariable, the mean stellar color in the reference field < H − K > ref can be used as the zero point of color excess. The color excess distribution toward globules can be derived by subtracting < H − K > ref from the observed H − K of stars. Since the extent of each analyzed field is relatively small (∼ 7 ′ × 7 ′ ), our assumption on invariant stellar population seems plausible. To measure color excess distribution, we arranged circular cells (∼30 ′′ diameter) spaced at 9 ′′ (the same alignment as used in the star count analysis) on each observed field. We calculated the mean stellar color in each cell and derived their color excess as where (H − K) i is the color index of the i-th star in a cell, and N is the number of stars falling in a cell. H − K ref is the mean color of stars in the reference field (see, Fig. 2). We converted the H −K color excess to A V using the reddening law of Rieke & Lebofsky (1985) as We then smoothed the resulting A V map in the same manner as used in the star count analysis.
-9 - The resulting A V maps are shown in Figure 2. Typical uncertainty of A V in the reference field (A V ∼ 0 mag) is estimated to be ∼0.6 mag.

Molecular Line Parameters
We derived the observed molecular line parameters, LSR velocity, line width, and optical depth, for each globule. We observed N 2 H + (J = 1 → 0) for FeSt 1-457 and C 18 O (J = 1 → 0) for the other globules. We show the spectrum observed at the center of each globule in Figure 3. For the N 2 H + spectrum fitting we took into account their hyperfine structure and line optical depth effect (using the same procedure as described in Tatematsu et al. 2004). The intrinsic relative intensity and rest frequency of hyperfine components are taken from Tiné et al. (2000) and Caselli et al. (1995), respectively.
The observed line width ∆V obs of the globule ranges from ∼ 0.2 to ∼ 0.5 km s −1 ( Table 3). The kinetic temperature T k of low mass dense cores is known to be ∼ 10 K (e.g., Benson & Myers 1989). Thus, it is natural to attribute the excess in the observed line width with respect to the thermal line width (0.13 km s −1 for C 18 O and N 2 H + ) to the contribution from non-thermal turbulence in globules. We estimated non-thermal line width ∆V NT of each globule by assuming that T k is equal to 10 K. ∆V NT can be obtained from ∆V 2 NT = ∆V 2 obs − 8ln2kT k /m obs , where m obs is the mass of the observed molecule. The thermal line width ∆V T = (8ln2kT k /m) 1/2 for T k = 10 K is 0.443 km s −1 , where m is the mean molecular weight (2.33 amu). Derived line parameters are listed in Table  3. We found that thermal line width for the mean molecular weight is larger than non-thermal line width for most of the observed globules. It is likely that observed globules are dominated by thermal support with smaller contribution from the non-thermal turbulence. The thermal and nonthermal components of the line width can be characterized by a temperature as T k = ∆V 2 T m/k8ln2 and T NT = ∆V 2 NT m/k8ln2, respectively. We derived T eff = T k + T NT for each globule in Table 3 as a cloud-supporting temperature which includes non-thermal contribution. In order to estimate the distance to the globules, we will use T eff in combination with the Bonnor-Ebert model fitting of globules in §4.2.2. In this paper, we only use the information on molecular line width. Detailed studies such as comparison between dust extinction and molecular distribution will be reported in a subsequent paper.

Bonnor-Ebert Sphere Model
In order to model the internal density structure of Bok globules, we used the Bonnor-Ebert sphere model (Bonnor 1956;Ebert 1955), which is a modified Lane-Emden equation describing a pressure-confined, self-gravitating isothermal gas sphere in hydrostatic equilibrium. The density structure of a Bonnor-Ebert sphere can be obtained by solving the following equation: -10where φ is the logarithmic density contrast, φ(ξ) = − ln(ρ/ρ c ); ρ and ρ c are the volume density and central volume density, respectively. ξ is the dimensionless radial parameter, ξ = (r/C s ) √ 4πGρ c ; r, G, and C s are the radius, gravitational constant, and isothermal sound speed, C s = (kT /m) 1/2 , where k, T , and m are the Boltzmann constant, kinetic temperature, and mean mass of the molecule (2.33 amu), respectively. Under the standard boundary conditions, φ(0) = 0 (i.e., ρ = ρ c at r = 0) and dφ(0)/dξ = 0 (i.e., dρ/dr = 0 at r = 0), equation (8) can be solved with numerical integration. If the gas sphere is confined by external pressure P ext at the core boundary radius R, the solutions of the Bonnor-Ebert sphere can be characterized by the dimensionless radius ξ max = ξ(r = R) as ξ max is a stability measure of the gas sphere against the gravitational collapse (Bonnor 1956;Ebert 1955). Critical state is achieved at ξ max = 6.5 corresponding to the density contrast of ρ c /ρ R = 14.
For a solution with ξ max > 6.5, the equilibrium state is unstable to the gravitational collapse. The Bonnor-Ebert sphere mass M and the external pressure P ext are calculated as where α = C s /(4πGρ c ) 1/2 .

Derivation of Column Density Profile
For a comparison of observational data (i.e., A V ) with a theoretical Bonnor-Ebert density distribution, radial A V profiles for each globule were constructed. We set annuli spaced at 9 ′′ around the center of the core that is the A V intensity peak (see Table 4), and calculated averaged A V values at each annulus (see Figure 4). For apparently elliptical globules (CB 87, CB 161, and CB 184) we performed ellipse fitting to the A V distribution, and "equivalent" radius, r = √ ab, was substituted for radius where a and b are the respective semi-minor and semi-major axes of the ellipse for the cores. Fitted ellipticity and position angle for the globules are listed in Table 4. It is possible that the observed globules with circular projected shapes are actually pole-on prolate ellipsoids. In this case, the shape of the observed A V profile should be steeper than that for the intrinsic shape. Since there is no effective method for estimating how ellongated a globule is along the line of sight, we assumed that the line of sight length of globules is the same as the projected diameter.
We manually masked extinction features likely to be unrelated to the globules (e.g., the other clouds on the same image, diffuse streaming feature; see Fig. 2), and masked regions were omitted in the circularly averaging procedure. Furthermore we recalibrated the relative zero point of extinction for the analyzed globules by subtracting background extinction (i.e., A V offset value). The A V offsets have been derived from averaging A V within the annulus of 150 ′′ to 200 ′′ radius, which is sufficiently larger than the primary extinction feature of globules. Actually, the subtracted values as A V offset are small (typically ∼ 0.5 mag), and it solely depends on where we set the reference fields in the star count or color excess analysis. In order to convert A V to the H 2 column density, we used the relationship between the optical color excess E B−V and the hydrogen column density, N (HI + H 2 )/E B−V = 5.8 × 10 21 cm −2 mag −1 (Bohlin, Savage, & Drake 1978), where N (HI + H 2 ) = N (HI) + 2N (H 2 ). This relationship is based on the UV absorption line measurements toward the stars of slightly obscured lines of sight (E B−V < ∼ 0.5 mag). By assuming that (1) all hydrogen in the cloud is in a molecular form and that (2) the ratio of total-to-selective extinction is equal to the standard ISM value of 3.1 (e.g., Whittet 1992), we obtained N (H 2 )/A V = 9.4 × 10 20 cm −2 mag −1 , and used it to derive the H 2 column density. The resulting column density profiles are shown in Figure 4.

Derivation of Model Parameters
In order to produce Bonnor-Ebert model density profiles, we performed numerical integration of equation (8). Since the Bonnor-Ebert sphere is characterized by a single parameter ξ max (Eq. [9]), best-fit ξ max for each globule can be determined by comparing model column density profiles with those obtained from observations. We integrated the model density profile along the lines of sight on grids to calculate the column density profile. In order to match the resolution between the model and the observed column density data, the model column density on the ∆α-∆δ plane was convolved with the beam used in the star count or color excess analysis (cylinder-shaped beam convolved with the Gaussian). In the χ 2 -procedure, the radial profile was fitted from the minimum radius (center) to the radius where the observed A V value reaches ∼ 1 mag (detection limit). In some cores (Lynds 495 and Lynds 498), deviation of observed profiles from the Bonnor-Ebert sphere is found at diffuse outskirts (A V < ∼ 1 mag), which was excluded in the fitting. The χ 2 -fitting results are not sensitive to the selection of the maximum radius for fitting part unless there is no significant kink in the column density profile.
In Figure 4, we show the results of the Bonnor-Ebert fit to the radial column density profile of the globules. The dots and error bars represent the average N (H 2 ) values at each annulus of 9 ′′ intervals and the rms deviation of data points in each annulus, respectively. The solid line denotes the radial column density profile of the best-fit Bonnor-Ebert sphere which is convolved with beam used in the A V measurements. The dashed line denotes the Bonnor-Ebert model profile before the convolution. In Table 4, we list the derived Bonnor-Ebert parameters for globules including previous reports on four sources (Barnard 68, Barnard 335, Coalsack Globule II, and Lynds 694-2). The error of each parameter in the table represents 1σ deviation in the χ 2 density profile fitting. We discuss possible systematic uncertainties in the derivation of ξ max in the Appendix.
Recently, Hotzel et al. (2002b) and Lai et al. (2003) presented detailed discussion on the relation between the physical parameters of the Bonnor-Ebert sphere. Substituting the relations of (9), we can obtain d −1 T =constant (Lai et al. 2003). The distance d or temperature T of globules should be accurately determined from observations. Since there were no accurate methods for measuring distance to the starless, compact, and nearby (d < ∼ 500 pc) globules except for the limited cases (e.g., interaction with HII region of known distance, existence of foreground stars toward them), the distances listed in Table 1 serve -12as the first reliable values. If we use the initially assumed distance d ini , the necessary temperature for the Bonnor-Ebert sphere T BE can be determined. If we assume that T eff (= T k + T NT ) from our molecular line observations (see §3.3) can approximate effective supporting temperature of the cloud ‡ , the corrected distance to the globule, d = (T eff /T BE )d ini , can be determined. Finally we derived all the physical parameters of best-fit Bonnor-Ebert spheres (e.g., radius, central density, and mass) for the above "distance-assumed" and "temperature-assumed" cases as listed in Table 5. We note that the "temperature-assumed" distance to the starless globule FeSt 1-457 is estimated to be ∼ 70 pc, which would make this globule one of the nearest dark clouds.

Discussion: Stability and Evolution of Globules
In this section, we discuss the internal structure and stability of globules by using the data of the ten globules from our observations and four from the literature.

Physical Properties of Globules
Though we have estimated physical properties of globules by assuming the Bonnor-Ebert model, there are actually other theoretical models (e.g., Shu 1977;Larson 1969;Penston 1969;McLaughlin & Pudritz 1997;Plummer 1911). Bacmann et al. (2000) presented that a Bonnor-Ebert-like model (a finite size sphere with inner uniform density region and ρ ∝ r −2 outer envelope) fits well with the column density profile of globules obtained from mid-infrared 7 µm absorption measurements, but that a singular-isothermal-sphere (SIS) or singular logotrope sphere cannot reproduce the observed profile well. By considering the resolution of our extinction maps (∼30 ′′ ) and possible geometric error of globules (e.g., deviation from spherical symmetry), it is not possible to make meaningful comparisons between the Bonnor-Ebert model and other Bonnor-Ebert-like models (see, e.g., Figure 2 in Harvey et al. 2003, which compares the Bonnor-Ebert model with a Plummer-like model). Here, we simply assume the Bonnor-Ebert model.
As described in the previous section, the dimensionless radial parameter ξ max determines the shape of a Bonnor-Ebert density profile as well as the stability of the equilibrium state against the gravitational collapse. We plot globules on the ξ max versus density contrast diagram in the upper panel of Figure 5. Since ξ max has one-to-one correspondence with density contrast, all the Bonnor-Ebert spheres are distributed along the solid curved line. We found that more than half of starless globules (7 out of 11 sources) are located near the critical state, ξ max = 6.5 ± 2. Thus, we suggest that a nearly critical Bonnor-Ebert sphere characterizes the typical density structure of starless globules. The remaining starless globules show clearly unstable states (ξ max > 10). When we divide starless globules into two groups with respect to the critical line (ξ max = 6.5), there are three stable starless globules and eight unstable starless globules. The majority of starless ‡ Theoretically, supporting force in a Bonnor-Ebert sphere is purely thermal pressure, but actual clouds also have non-thermal support, e.g., turbulence and magnetic fields, which are widely observed in molecular clouds (see, e.g., Larson 2003). We assume that non-thermal support can be included in the model by replacing the isothermal sound speed Cs = (kT k /m) 1/2 in equation (9) by the effective sound speed C eff = (kT eff /m) 1/2 . globules is located in the unstable states (ξ max > 6.5) if the uncertainties in ξ max are not taken into account. We also found that all the star-forming globules have larger ξ max values (> 10) than that of the critical Bonnor-Ebert equilibrium state, which is consistent with the fact that they have started gravitational collapse. In the lower panel of Figure 5, we show the fitting error in ξ max for each globule. Filled gray circles denote the fitting results for some globules with considerable features of ambient extinction if the masked region on the A V map (see Fig 2) is included in the derivation of the column density profile. In addition to the ξ max diagram in Figure 5, we made histograms of logarithmic density contrast for starless and star-forming globules ( Figure 6). Similar characteristics as descrived above can also be seen in the histograms.
To see the relationships between derived physical properties of globules in Tables 4 and 5, we made correlation diagrams between density contrast and the other physical parameters of globules in Figure 7. We consider that the horizontal axis (density contrast) represents the evolutionary states of globules, because globules should evolve toward higher central condensation to form stars. We can see from panels (e) and (f) in the figure that temperature T eff and external pressure P ext of globules appear roughly constant regardless of the density contrast. If we assume constant T eff and P ext , the Bonnor-Ebert model parameters for each density contrast value can be constrained. The broken line in each panel is the relationship for the Bonnor-Ebert spheres with a constant T eff of 13.2 K and P ext of 5.3 × 10 4 K cm −3 , which are the median values for the globule samples. The value of P ext is larger than the previously reported ISM pressure in the solar vicinity of ∼ 1.8 × 10 4 K cm −3 (McKee 1999).
The overall trend in the data point distributions in Figure 7 can be approximated by Bonnor-Ebert spheres with a constant "effective" temperature and a constant external pressure. The relations for the radius and mass (panels [c] and [d]) remain roughly constant against the density contrast, suggesting that these quantities are not sensitive to the evolutionary states of globules. This result agrees with the N 2 H + observations of molecular cloud cores in Taurus (Tatematsu et al. 2004), in which there is no remarkable difference in core radius and mass between starless cores and cores with stars. Since the density at core boundary is derived from P ext /T eff , the density at core center n c is directly obtained for each density contrast value and shows a proportional relationship against the density contrast as shown in panel (b). Accordingly, the density contrast and A V show positive correlation. The data points in panel (a) are located slightly lower than the broken line because of a beam dilution effect in the A V measurements with ∼ 30 ′′ resolution. For molecular cloud cores located in regions of higher pressure (e.g., dark cloud complex) than in the solar vicinity, the Bonnor-Ebert relationship in panels (a) and (b) should be shifted upward corresponding to the larger value of P ext /T eff . We can see from panel (a) that the density contrast value is simply estimated from the A V (column density) measurements under the assumption of the Bonnor-Ebert sphere with a constant T eff and P ext .

On the Selection Bias of the Globule Samples
Here we briefly consider the selection bias in our globule samples. The peak A V values for our observed globules range from ∼ 6.5 to 41 mag. We selected globules that were opaque at optical wavelengths, i.e., with few stars in the background of the cloud on the Digitized Sky Survey image.
-14 -These sources are expected to have A V greater than at least 5 mag, and expected to be sites of future star formation, i.e., so-called "dense cores". The lower limit of the peak A V for the globules may be determined by the criterion in our target selection. Thus, it is possible that the diffuse cores of low column densities are not included in our globule sample as in the case of previous studies on dense cores (e.g., Visser et al. 2002;Jijina et al. 1999). We note that the lower limit A V of ∼ 5 mag corresponds to a low density contrast of ∼ 3 if we assume a constant T eff and P ext used in the previous section. The density contrast value is much less than the critical value of 14.

Implications for the Stability of Globules
It is interesting that there are starless globules showing larger density contrast than that of the maximum value (∼14) for stable Bonnor-Ebert spheres (see § 5.1), because unstable equilibrium states should not last long without any extra stabilizing force such as magnetic pressure and/or turbulent pressure (e.g., McKee & Holliman 1999). Since unstable equilibrium states cannot last long, the globules of large density contrast should not exactly be the Bonnor-Ebert sphere even if they mimic the Bonnor-Ebert density profile. Similar results have also been reported from the mid-infrared extinction study of starless dense cores (density contrasts are generally ∼10−80; see Table 3 of Bacmann et al. 2000). If we consider a simple case in which the non-thermal pressure components have the same radial dependence as the thermal pressure, the theoretical Bonnor-Ebert model can be modified by replacing its kinetic temperature term with T eff as described in § 4.2.2. In this case, the shape of the Bonnor-Ebert density profile as well as the stability criteria (i.e., the value of ξ max or maximum density contrast) hold unchanged compared with the purely thermally supported case. This modification cannot provide a realistic explanation for the stability of the large density contrast globules, suggesting either that they are already collapsing toward higher central condensation as discussed in § 5.5, or that the stabilizing mechanisms accounting for their large density contrast are more complex. The Bonnor-Ebert assumption for the globules is not valid in these cases. We note that our result is in marked contrast to one in which most of the clumps observed at 850 µm dust emission in the Orion B molecular cloud complex are in stable equilibrium based on the Bonnor-Ebert analysis (Johnstone et al. 2001).
A non-isentropic multi-pressure polytrope model (McKee & Holliman 1999) may provide successful interpretation allowing larger density contrast than that of the critical Bonnor-Ebert sphere. Curry & McKee (2000) developed a composite polytrope model which also allows large density contrast. These two models provide a stable equilibrium solution of arbitrarily large density contrast by adjusting the contributions from the non-thermal pressure components, i.e., turbulent and magnetic pressure. Galli et al. (2002) studied the effects of the interstellar radiation field (ISRF) on the cloud stability. They reported that a cloud's maximum density contrast increases with increasing intensity of the ISRF if the stabilizing effect of external heating of the cloud as well as produced temperature gradient are considered. These models are possible alternatives to the Bonnor-Ebert sphere providing larger value of maximum density contrast (> 14). The shape of their density profiles have characteristics similar to those of the Bonnor-Ebert sphere (see, Fig. 3 and 4 in Galli et al. 2002, Fig. 10 in McKee & Holliman 1999), that is a flat inner region surrounded by steeper envelope decreasing with power-law, in particular ρ ∝ r −α where α ∼ 2. It seems difficult to distinguish them only from the observed shape of the density profile. In other words, each model should provide similar density profiles and corresponding density contrast values for the same observational dataset. Detailed observations of velocity and temperature structure of dense cores are needed to determine the most realistic model.
Major differences among the models discussed above are the supporting mechanisms of the cloud and the stability criterion, i.e. the value of the critical density contrast. The actual density contrast value for globules can be obtained through the Bonnor-Ebert fit (see Table 4), but its critical value is not always known (i.e., model dependent). Initial density structure of globules right before the collapse can be determined if we constrain the value of its critical density contrast. Initial density conditions for star formation will provide a good starting point for the theoretical calculations. For example, given initial density profile of a globule, time evolution of mass accretion rate can be calculated (e.g. Whitworth & Ward-Thompson, 2001), which is an important parameter for constraining the mass of a newly formed star.
In order to estimate the value of critical density contrast from observations, it should be worth investigating the frequency distribution of density contrast for a number of starless globules. Since unstable equilibrium states should not be long sustainable, we expect that the detection probability for stable globules is much higher than that for unstable ones. Thus, it is likely that the frequency distribution of density contrast shows a sharp decrease toward larger density contrasts at a certain bin, which should correspond to the density contrast critical for collapse. Though the number of our samples (11 sources) is not large, we present the histogram of density contrast for the starless globules in Figure 6. The density contrast of the starless globules is mostly populated near the critical Bonnor-Ebert sphere and decreases toward the larger density contrast, suggesting that the actual critical density contrast of globules is not very different from that for the Bonnor-Ebert sphere. A statistical study with a larger number of samples is our future plan.

Starless Globules with Stable Bonnor-Ebert Solutions
In Table 4, three out of fourteen globules are fitted as a stable equilibrium state of the Bonnor-Ebert sphere. If they are truly stable, how can stable globules form stars? Stable globules might increase their density contrast quasi-statically toward the onset of dynamical collapse with gradually increasing external pressure (e.g., Hennebelle et al. 2003) and/or with decreasing internal turbulent pressure (Nakano 1998). For the case of rapid increase of external pressure, core collapse triggered by a compression wave should be taken into account (e.g., Hennebelle et al. 2003Hennebelle et al. , 2004Motoyama et al. 2004). Our observations of globules, however, show neither systematic decrease in T eff nor increase in P ext (panels e and f of Fig. 7) toward larger density contrast. It seems that non-thermal turbulence is already mostly dissipated even in the stable globules. Though there is only small room for the dissipation of turbulence, a slight decrease of turbulence may initiate the collapse of stable globules as suggested by Nakano (1998). The decrease of T eff (i.e., line widths) necessary for stable globules to start gravitational collapse is small. This is consistent with the result that there is no apparent tendensity of systematic decrease in T eff against density contrast (panel e of Fig.  7). In the star forming regions where non-thermal motions are more prominent than globules, e.g., -16 -GMCs, the dissipation of turbulence in dense cores is reported through the virial analysis (e.g., OMC-2/3 region in Orion: Aso et al. 2000).
Since possible evolutionary processes of stable globules are expected not to be dynamical, their lifetime should be long compared with unstable ones. However, in our starless globule samples, the fraction of stable globules (density contrast is less than 14) is small (3 out of 11 sources). The observed starless globules are mostly populated near the critical Bonnor-Ebert sphere, and decreases toward the smaller density contrast (Figure 6). This result is inconsistent with an idea that all of nearly critical globules originated from stable low-density-contrast globules; if the idea is the case, the number of the observed low-density-contrast globules should be larger than or comparable with that of the nearly critical ones. Thus, it is most likely that most of the nearly critical globules which we observed are not evolved from stable globules. The origin of nearly critical globules can be related to the formation mechanism of globules. Further studies on this subject are desirable.

Slow Collapse of A Nearly Critical Bonnor-Ebert Sphere
After gravitational collapse begins, globules evolve toward higher central condensation. If we assume an extreme case in which all the globules fitted as unstable Bonnor-Ebert equilibrium states are already collapsing, the difference in ξ max is attributed to different evolutionary stages of globules. Since the Bonnor-Ebert spheres provide a good fit to starless globules, it is meaningful to consider the gravitational collapse of the Bonnor-Ebert sphere (e.g., Foster & Chevalier 1993;Ogino, Tomisaka, & Nakamura 1999). We expect that the critical Bonnor-Ebert sphere can approximate the initial condition of collapse, because a large fraction of the starless globules are located near the critical state, ξ max = 6.5 ± 2, as shown in § 5.1.
In order to make quantitative comparisons of the density structures between our observations and a collapsing gas sphere starting from a nearly critical Bonnor-Ebert sphere, we used the result of theoretical calculations from Aikawa et al. (2005). In their calculation, the central density n c (H) and kinetic temperature T of the critical Bonnor-Ebert sphere were set as n c (H) = 2 × 10 4 cm −3 and T = 10 K, respectively. No cloud rotation, magnetic fields, or non-thermal turbulence were considered, and a fixed boundary condition of zero velocity was adopted at the outermost radius of R = 0.2 pc, which is slightly larger than the critical radius of an equilibrium sphere (see, Aikawa et al. 2005). We note that the initial condition parameters of n c (H), T , and R, are consistent with the observed physical properties of globules around ξ max = 6.5 (see , Table 4 and 5). The initial density distribution was multiplied by a constant factor α of 1.1 in order to initiate collapse. The slight density enhancement of α = 1.1 is preferable to approximate the collapse from a nearly critical Bonnor-Ebert sphere. The time scale of the collapse is several times the free-fall time, ∼ 1.2 × 10 6 yr (∼ 4 × t ff ), which is comparable with the observational lifetime of starless dense cores (∼ 10 6 yr; e.g., Visser et al. 2002;. Figure 8 shows the evolution of a radial column density profile of the collapsing sphere for the case of α = 1.1 described above. The gray solid line denotes initial column density profile (t = 0 yr), and the black solid, dotted, and broken lines denote density profile at specific times of 4.4×10 5 , 8.9 × 10 5 , and 1.1 × 10 6 yr, respectively. We found that the column density profiles of the collapsing -17sphere are well fitted with the static Bonnor-Ebert spheres (equilibrium solutions) in the range of ξ max < 25 within the relative difference of 2%. The plus symbols show the best-fitting result for each collapsing profile with the Bonnor-Ebert sphere, which are labeled with corresponding ξ max values; for example, the density profile of the collapsing sphere at 4.4 × 10 5 yr mimics the unstable Bonnor-Ebert sphere of ξ max = 7.6. This result suggests that observed starless globules with large density contrasts are gravitationally collapsing objects. Figure 9 shows the relationship between ξ max versus density contrast. The solid line is for the Bonnor-Ebert sphere, and the dotted line is for the collapsing sphere which plots the best-fitting ξ max values against the density contrast of the collapsing sphere. The plus symbols denote the elapsed time of the collapsing sphere at 10% intervals of the total collapse time. Since the free-fall time scale increases with decreasing density, a collapsing gas sphere starting from a nearly critical Bonnor-Ebert sphere should mimic marginally unstable Bonnor-Ebert spheres for a long time. Thus, this is qualitatively consistent with the fact that the ξ max for the starless globules is mostly populated near the critical state (7 out of 11 sources) and the remaining starless globules show clearly unstable states (ξ max > 10; 4 out of 11 sources) as shown in § 5.1. From the model prediction, we note that 50% of collapsing globules should be located in ξ max = 6.5 − 8.5, and the rest of the globules should have ξ max larger than 8.5. The fraction is consistent with our result that half of the starless globules with ξ max > 6.5 are located in ξ max = 6.5 − 8.5 (4 out of 8 sources, see Fig. 5). In Figure 10, we demonstrate the histogram of logarithmic density contrast for starless globules (the same as Figure 6) and the time evolution of the model collapsing sphere (solid line with dots). Each dot represents the residence time of a collapsing sphere measured at 0.2 intervals in the logarithmic density contrast. The histogram of the observed globule in the unstable region is consistent with that from model calculations of a collapsing sphere.
As discussed above, the density structure of the observed globules can be explained in terms of the slow collapse of a nearly critical Bonnor-Ebert sphere. The anticipated inward velocity from the model calculation remains subsonic § for a long time; the peak inward velocity is less than ∼ 0.1 km s −1 for ∼ 10 6 yr, which is ∼ 80% of the total collapse time (see, Fig. 1 of Aikawa et al. 2005). Previous reports on the observations of CS (J = 2 → 1) (blue) asymmetry profile, indicative of gravitational collapse in dense cores, show that infall candidates constitute 25 % of the starless core samples (Lee, Myers, & Tafalla 1999). It is possible that the rest of the starless cores have a slow collapse motion which is hard to detect. Radio molecular line observations of infall asymmetry profiles with sufficiently high velocity resolution (e.g. ∼ 0.01 km s −1 ) are necessary to confirm the scenario of slowly collapsing globules. Among the starless globules in Figure 5, the gas inward velocity for Lynds 694-2 (ξ max = 25 ± 3) was measured as 0.05 − 0.07 km s −1 (Lee, Myers, & Tafalla 2001;Lee, Myers & Plume 2004), which is consistent with the prediction from the collapse model (α = 1.1).

Summary
We have carried out a near-infrared (J, H, and K s ) imaging survey of ten Bok globules using the infrared camera SIRIUS on the IRSF 1.4 m telescope in South Africa. A V distributions for the ten globules were derived through measurements of the H − K color excess or the stellar density at J, H, and K s (star count). The radial column density profile for each globule was analyzed using the Bonnor-Ebert sphere model, which well fits the profile within observational uncertainties, and physical properties of globules, e.g., size, central density, temperature, mass, external pressure, and center-to-edge density contrast, were derived. The dimensionless radial parameter ξ max determines the shape of a Bonnor-Ebert density profile as well as the stability of the equilibrium state, e.g., the solution of ξ max > 6.5 is unstable to the gravitational collapse. We investigated the stability of globules on the basis of ξ max for ten globules from our observations and four globules in the literature. In addition to the near-infrared imaging, we have carried out radio molecular line observations toward the same ten Bok globules using the 25-element focal-plane SIS receiver BEARS on the 45 m telescope of Nobeyama Radio Observatory (NRO). We measured line width of the globules for independent measurements of the "effective" temperature including the contribution of turbulence. It was confirmed that most of the globules are dominated by thermal support. Since the distance-dependent parameters of the Bonnor-Ebert sphere were derived using the initially assumed distance to the globules, derived best-fitting temperature was affected by the distance assumption. We constrained the distance to the globules by comparing the temperature from the molecular line width with the temperature from the Bonnor-Ebert fitting, and rescaled distance-dependent physical quantities of globules.
(1) We found that more than half of the starless globules (7 out of 11 sources) are located near the critical state, ξ max = 6.5 ± 2. Thus, we suggest that a nearly critical Bonnor-Ebert sphere characterizes the typical density structure of starless globules, and it approximates the initial condition of gravitational collapse.
(2) We found that four out of eleven starless globules show clearly unstable states (ξ max > 10). Since unstable equilibrium states should not be long sustained, we expect that they are already collapsing toward higher central condensation or that extra force (e.g., magnetic and/or turbulent pressure) accounting for large ξ max stabilizes the globules. It was also found that all three starforming globules have unstable solutions of ξ max > 10, which is consistent with the fact that they have started gravitational collapse.
(3) We investigated the collapse of the Bonnor-Ebert sphere from a nearly critical state using the model calculation of Aikawa et al. (2005), and found that the column density profiles of the collapsing sphere mimic those of static Bonnor-Ebert spheres (unstable equilibrium solutions). By relating ξ max to the collapsing sphere at a specific time, the evolutionary state of globules can be interpreted, and the detection probability of each ξ max value can be predicted from the model calculation. Since the evolutionary timescale decreases with increasing density, the collapsing sphere resembles a marginally unstable Bonnor-Ebert sphere for a long time. It was found that the frequency distribution of ξ max for the observed starless globules is consistent with that from model calculations of the collapsing sphere.    Table 1), Second row: value from molecular line observations (T eff ; see, Table 3). b First row: same value as listed in Table 1, Second row: calibrated value by using T eff to scale the fitted Bonnor-Ebert model parameters (see, §4.2.2). c The values in the first and second rows are derived using the distances listed in the third column.    The dots and error bars represent the average N (H 2 ) values at each annulus arranged at intervals of 9 ′′ and the rms dispersion of data points in each annulus, respectively. The solid line denotes the best-fit Bonnor-Ebert profile which is convolved to the observing resolution (typically ∼30 ′′ ) with the same beam as used in the A V measurements. The dashed line denotes the Bonnor-Ebert profile before the convolution. The gray dot-dashed line denotes the 1 σ deviation of the column density measured in the reference field. .-The solid line shows the relation between ξ max and density contrast (center-to-edge density ratio). The location of our 10 globules and previously reported 4 globules (Barnard 335, Harvey et al. 2001;Barnard 68, Alves et al. 2001; Coalsack Globule II, Racca et al. 2002;Lada et al. 2004;Lynds 694-2, Harvey et al. 2003) are shown in the figure. Filled and open black circles denote starless and star-forming globules, respectively. Filled gray circles denote the fitting results when the masked region on the A V map (see Fig. 2) is included in the derivation of the column density profile. The vertical dashed line denotes critical state of ξ max = 6.5. In the lower panel of the figure, horizontal error bars denote 1σ error in the χ 2 density profile fitting. -Correlations between density contrast and the other physical parameters. The panels of (a) to (f) show the correlation diagrams for A V , central density, radius, mass, temperature, and external pressure, respectively. Filled and open circles denote starless and star-forming globules, respectively. The vertical gray line denotes the density contrast of the critical Bonnor-Ebert sphere. The plotted data points are the physical quantities after distance-correction (see, §4.2.2, Table 5). The physical quantities without distance-correction are shown as the symbol "−" in each panel. The broken lines represent the relationships for the Bonnor-Ebert spheres with a constant temperature T eff of 13.2 K and P ext of 5.3 × 10 4 K cm −3 .
-38 - Fig. 8.-Radial column density distribution of a collapsing gas sphere at some specific times (black solid, dotted, and broken line) started from the initial condition of a nearly critical Bonnor-Ebert density structure (gray solid line). Dotted and broken lines are shifted by 5 and 10 (10 21 cm −2 ) in the figure, respectively. Best-fit solution of the Bonnor-Ebert sphere for each line is also plotted as plus symbol.
-39 - Fig. 9.-The solid line shows the relation between ξ max and density contrast for the Bonnor-Ebert sphere (same as Fig. 5). The dotted line denotes the similar relation for a collapsing sphere starting from a nearly critical Bonnor-Ebert sphere (α = 1.1 model; Aikawa et al. 2004); the column density profiles of the collapsing sphere were fitted with static Bonnor-Ebert sphere solutions, and the bestfit ξ max values are plotted against the density contrast values of the collapsing sphere. The plus symbols denote the elapsed times from the onset of collapse at 10% intervals of the total collapse time (t col ).
-40 - Fig. 10.-The histogram of the logarithmic density contrast for starless globules (the same as Figure 6). The solid line with dots shows the time evolution of the model collapsing sphere; each dot represents the residence time of a collapsing sphere measured at 0.2 intervals in the logarithmic density contrast. The vertical axis on the right side shows the percentage of the residence time with respect to the total collapse time.