© ESO, 2012

1. Introduction

It is well known that for most globular clusters the relevant relaxation times are shorter than their age. Therefore, since the two-body relaxation processes have had time to act, these systems may be considered to be close to thermodynamical relaxation, with a distribution function that should be close to a Maxwellian.

The spherical King models (King 1966) are defined by a distribution function of this type, with a truncation in energy introduced to mimic the role of tides (for the relevant definition and a summary of the main properties, see Appendix A.1). In general, they are thought to be a physically justified simple tool to represent quasi-relaxed stellar systems, in particular well suited to describe globular clusters. In this context, King models are usually considered as the correct zeroth-order dynamical reference model. The available catalogs (e.g., see Harris 2010) record a number of structural properties of globular clusters that are derived from an interpretation of observational data based on the use of King models. The family of King models is also frequently used in dynamical studies aimed at identifying given physical mechanisms and represents a most popular choice of initial conditions for numerical simulations designed to investigate the dynamical evolution of star clusters (e.g., see Chernoff & Weinberg 1990; and, more recently, Gieles & Baumgardt 2008).

In reality, this simple physical picture suffers from a number of limitations, which become more evident now that much improved observations allow us to go well beyond a zeroth-order picture of globular clusters. Of course, numerical simulations can now model in a realistic way the important effects of discreteness, mass segregation, binaries, rotation, core collapse, and time-dependent tides (e.g., see Küpper et al. 2010), which are known to govern the evolution of globular clusters. Still, in relation to the general issue of relaxation it would be interesting to test further to what extent King models (or alternative idealized models) actually work as a useful paradigm to represent small, quasi-relaxed stellar systems.

From the observational point of view, there remains indeed one unsatisfactory point, which this paper would like to address, at least on the basis of the currently available data. This is that the application (and thus the test) of the spherical King models is frequently carried out, in practice, only in relation to the available photometric profiles (see Trager et al. 1995; McLaughlin & van der Marel 2005), without the corresponding tests on the associated kinematical profiles, in contrast to what is routinely done in the studies of early-type galaxies. Some studies of Galactic globular clusters, constrained simultaneously by density and velocity dispersion profiles, are actually available, mainly based on multi-mass Michie-King models (e.g., see Gunn & Griffin 1979; Meylan et al. 1995; Meylan & Mayor 1991, for M3, ω Centauri, and NGC 6397, respectively), but a systematic and homogeneous investigation is still missing. Note that, on the galactic side, deep investigations of these issues, starting with the mid 1970s, have led to the remarkable discovery that bright ellipticals are generally supported by anisotropic pressure and contain significant amounts of dark matter inside the effective radius R_e (e.g., see Chap. 24 in Bertin 2000). From the study of elliptical galaxies it has also been learned that structurally different models (as diagnosed by their kinematics or characterized by their virial coefficients) may have remarkably similar photometric profiles (e.g., see Appendix D of Bertin et al. 2002).

In this paper, we wish to check whether indeed King models perform better in more relaxed systems, by studying the combined photometric and kinematical profiles for a sample of globular clusters. We will refer to the surface brightness profiles collected by Trager et al. (1995) (supplemented by more recent data, when available), in order to deal with a homogeneous sample, and limit our discussion to one-component dynamical models (leaving aside the issues of core collapse and mass segregation; we will thus exclude from our sample clusters with evidence of core collapse). Our sample will then be basically defined by the requirement that a sufficient number of kinematical data-points is available from the literature, so as to define a reasonably accurate and radially extended kinematic profile for a test of a dynamical model at the global level.

Many large globular clusters have very long relaxation times. Therefore, we have decided to model the same data also by means of models (the f^(ν) models; see Bertin & Trenti 2003) explicitly constructed for the context of violently relaxed elliptical galaxies. In other words, we wish to check whether less relaxed clusters tend to conform to the picture of formation via incomplete violent relaxation, which has the characteristic signature of radially-biased pressure anisotropy for less bound stars. The use of the f^(ν) models is preferred to other options (for example to the use of King-Michie models, Michie 1963), because these models are based on a detailed physical justification and have been shown to perform well both in relation to the observations of bright ellipticals and to the properties of the products of incomplete violent relaxation found in numerical simulations of collisionless collapse (over a range of nine orders of magnitude in the computed density profiles, with an excellent fit to the properties of the generated pressure anisotropy profiles; see Trenti et al. 2005). For a fair comparison with the King models we should have referred to truncated f^(ν) models; for the present simple exploratory investigation, and in order to keep the comparison between models characterized by the same number of parameters, we have decided to use the non-truncated f^(ν) models.

The paper is organized as follows. In Sect. 2 we introduce the sample of globular clusters selected for this study. In Sect. 3 we describe the available data sets and the procedure we followed to obtain the radial profiles used in our analysis. In Sect. 4 we show the results of our work and in Sect. 5 we draw our conclusions. A brief description of the characteristics of the dynamical models is provided in Appendix A; in Appendix B we describe the adopted fitting procedure.

2. The selected sample

In this paper we wish to analyze globular clusters characterized by different relaxation conditions, measured by the central (core) relaxation time T_c as listed in the Harris catalog (Harris 2010). We decided to order globular clusters according to the core relaxation time, rather than to the half-mass relaxation time, because in general T_c is less model-dependent (as pointed out in Sect. 4.2). We checked that, by ordering the sample with respect to the half-mass relaxation time, some changes in the composition of the three classes identified below would occur (this fact can be easily seen by inspecting the values of these time scales for the globular clusters listed in Table 1), but with minor effects in relation to the general conclusions of the paper.

Looking at the values of the relevant relaxation times recorded in the available catalog for the globular cluster system of our Galaxy, we see that the distinction between relaxed and partially relaxed globular clusters is not sharp. Therefore, we introduce a simple criterion to classify globular clusters according to their relaxation state. We consider three relaxation classes: relaxed globular clusters, for which log T_c < 8 (first class); globular clusters in an intermediate relaxation condition, for which 8 < log T_c < 9 (second class); partially relaxed globular clusters, for which log T_c > 9 (third class; T_c is expressed in years).

The sample of globular clusters has been selected on the basis of the following criteria: (i) we exclude post-core-collapse globular clusters, that is, we reject the clusters labeled as post-core-collapse by Harris (2010). The reason is that we wish to test King models on the global scale, avoiding the subtle modeling issues that characterize the central regions of these systems, especially if phenomena a priori known to go beyond the King modeling are involved. (ii) We choose globular clusters for which an accurate and extended surface brightness profile is available. (iii) We select clusters for which at least 140 stellar radial (line-of-sight) velocities have already been measured. We impose a lower limit to the number of measured velocities, because we wish to extract from such data a reasonably well-defined velocity dispersion profile; the value of this limit is fixed in a way that allows us to include in the analysis globular clusters that belong to the different relaxation classes defined previously. (iv) We exclude from our list all the clusters that have less than 35 stellar radial velocities inside the projected half-light radius R_e (also called effective radius; the values of this quantity are reported in the Harris catalog, where the notation is r_h). We introduce this further requirement, because we wish to analyze velocity dispersion profiles that characterize the stellar systems on the largest radial extent. For bright elliptical galaxies the kinematical data-points inside R_e are usually the easiest to get, and often turn out to discriminate among very different models. To test how well the King models perform, the central regions are a natural ground for comparison with other models.

The most restrictive elements in identifying a significant sample of globular clusters are the requirements on the radial velocity data, because for only few globular clusters the desired data are available. Indeed, of the 28 Galactic globular clusters with a reasonable number of radial velocities (i.e., at least 40 line-of-sight velocity measures on the entire spatial extent of the cluster), 3 are flagged as post-core-collapse, 9 have less than 140 velocity data, and 3 have less than 35 data inside their projected half-light radius. In the Harris catalog, NGC 362 and NGC 7078 are indicated as possible post-core-collapse clusters, but we decided to keep them because, according to their concentration parameter in the Harris catalog, it is still possible to obtain an acceptable fit with King models. In this way, we are left with 13 globular clusters that match our selection criteria. In relation to the relaxation classes defined above, our set of globular clusters contains 5 well-relaxed clusters, 6 clusters in an intermediate relaxation condition, and 2 partially relaxed clusters.

To better characterize our sample in terms of the radial extent of their radial velocity data, we consider the ratio of the radius of the last kinematical point to the projected half-light radius¹, R_K/R_e, and the ratio of the radius of the last kinematical point to the truncation radius², R_K/r_tr. We judge the following values of the two ratios, R_K/R_e ≥ 3 and R_K/r_tr ≥ 0.8, to be satisfactory. All the selected globular clusters satisfy the first relation, and all but four globular clusters satisfy the second condition.

Table 1 gives the sample of selected globular clusters, listed in order of increasing core relaxation time log T_c. The first part of the table contains relaxed globular clusters, the second part those in an intermediate relaxation condition, and partially relaxed clusters are shown in the last part. For each object, the values of the adopted cluster distance from the Sun d_⊙ (expressed in kpc), the concentration parameter C, the logarithm of the core relaxation time log T_c, the logarithm of the half-mass relaxation time log T_M (where T_c and T_M, in other papers often indicated with the symbols t_rc and t_rh, are expressed in years) and the ellipticity e are recorded (as listed in the Harris 2010, catalog). In addition, the maximum rotational velocity v_rot (in km s^-1; the references for these values are Lane et al. 2011, for the globular clusters for which we use the kinematic data published in this paper, and Meylan & Heggie 1997, for the others), the number of velocity data-points available N_v, the number of velocity data-points inside the projected half-light radius N_e, the ratios of the radius of the outermost velocity point R_K to the projected half-light radius R_e (in the Harris catalog and in other papers often indicated as r_h) and to the truncation radius r_tr, and the sources of the kinematic data are given in the last columns.

For completeness, in Table 2 we list globular clusters with a fairly rich kinematical information available (line-of-sight velocity measurements for at least 40 stars), but rejected in the present analysis because they do not satisfy our selection criteria. The table is organized in three parts, thus showing which criterion is not met: in the first part we list post-core-collapse globular clusters, in the second part those with an insufficient total number of kinematic data (N_v), and in the third part those with an insufficient number of data inside the half-light radius (N_e). In each part of the table, globular clusters are listed in order of increasing number of total stellar velocity data N_v. Column entries are the same as in Table 1 (except for cluster distance, ellipticity, and rotational velocity, which are not recorded).

3. The data sets

3.1. The surface brightness profiles

To deal with a homogeneous sample, we decided to use the surface brightness profiles provided by Trager et al. (1995) (this is the same starting point of McLaughlin & van der Marel 2005). This choice guarantees that the profiles have been constructed with the same method, even though the actual data come from different sources.

For each globular cluster, the profile is composed of N_p photometric data-points, given by log R_i, the logarithm of the radius R_i, measured in arcsec, and by m_V(R_i), the V band surface brightness measured in mag arcsec^-2 at the radial position R_i. The data also include m_V,C(R_i), the surface brightness calculated with the Chebyshev polynomials, which provides an accurate approximation to the overall profile; m_V(R_i) − m_V,C(R_i), the Chebyshev residual; w_i, a weight that the authors assign to each measurement. All the surface brightness profiles considered in this paper are relative to the V band.

These data have to be properly corrected and treated before a comparison can be made with the theoretical models. First, we introduced an extinction correction, under the assumption that such extinction can be considered to be constant over the entire extent of each globular cluster. We calculated the extinction A_V from the reddening listed in the Harris catalog and obtained m(R_i) = m_V(R_i) − A_V, for i = 1,...,N_p. This is the only correction applied to the data: we assume that Trager et al. (1995) already removed any foreground and background contamination that could affect the measurements. Then, we followed the procedure described by McLaughlin & van der Marel (2005) to estimate the uncertainties δm_i on the data: starting from the weights, we calculated³δm_i = σ_μ/w_i (note that in the fits we use only the points with weights w_i ≥ 0.15 in the original profile, as suggested by McLaughlin & van der Marel 2005).

Usually, the less reliable parts in the profiles published by Trager et al. (1995) are the central regions. Therefore, we decided to combine these profile with the more recent and accurate surface brightness profiles by Noyola & Gebhardt (2006), when available. For the globular clusters NGC 104 and NGC 6341, we simply combined the profiles from the two sources; for NGC 6254 and NGC 7078, we decided to combine the two data sets by removing the points from Trager et al. (1995), when they do not match the profile by Noyola & Gebhardt (2006) (the removed points are inside 3.5″ and 14.5″, respectively). In the case of NGC 7078, it should be emphasized that, with this treatment, the profile changes significantly in the central regions, and the central slope becomes steeper.

In the case of NGC 5139 (ω Cen) we decided to add to the surface brightness profile the inner points that Eva Noyola kindly provided us (Noyola et al. 2008); still the number of data-points in our final composite surface brightness profile of this cluster is significantly smaller than that of the other clusters.

3.2. The velocity dispersion profiles

We divided the data (see Table 1 for detailed references) in several radial bins containing an equal number of stars; we chose the binning that represents the best compromise between having a rich profile and having accurate points (by increasing the number of bins, that is by decreasing the number of points per bin, the errors on the velocity dispersion increase). In principle, we might consider using unbinned data, to avoid loss of information (for the study of dwarf spheroidals, see Wilkinson et al. 2002). Here we preferred to follow the more traditional approach of constructing the associated one-dimensional profiles, a method that can be applied in a similar way to both kinematical and photometric data and that allows us to follow well-established fitting procedures used in the past (especially in studies of elliptical galaxies). For the clusters with less numerous kinematical data-points, we did experiment with changing the radial binning and checked the related consequences.

The method used to calculate the mean velocity and the velocity dispersion from stellar radial velocities is basically the one described by Pryor & Meylan (1993). We started by calculating the mean v_r and the dispersion σ_c for the entire set of velocities. The mean velocity represents the overall velocity of the entire cluster. Then, taking this value as a constant for the entire cluster, we calculated the line-of-sight velocity dispersion σ(R_i) and the related accuracy δσ_i inside the bins in which the data have been divided. For each bin, we indicate the distance from the center R_i as the mean of the radial positions of the stars that it contains. In this paper we ignore the possible presence of rotation and therefore consider the various kinematical data-points in each bin, after subtraction of the systemic velocity, to contribute only to velocity dispersion (random motions).

For the majority of the clusters in our sample, only one data-set of stellar radial velocities is available; in the following we discuss in detail the cases in which a composition of different data-sets has been performed or which require some additional comments. For NGC 104, two data-sets of radial velocities are available; we noticed that the data from Gebhardt et al. (1995) are more centrally concentrated than those from Lane et al. (2011)⁴. In order to have a complete sampling on the entire radial extent of the cluster, we decided to define a mixed data-set, composed of 499 data from Gebhardt et al. (1995), located inside 100″, and 2139 data from Lane et al. (2011), located outside that radius. The second case is cluster NGC 288, studied by Carretta et al. (2009) and by Lane et al. (2011). Since these papers publish the coordinates of each star, we were able to single out the stars in common between the two data-sets: for the stars in the overlap, velocity measures by Carretta et al. (2009), being more accurate, have been preferred. Finally, we excluded three stars for which the value of the velocity deviates by more than 4σ_c from the mean radial velocity v_r, obtaining a final sample of 171 data.

For the globular clusters NGC 6121, NGC 6656 and NGC 6809, in addition to recent data from Lane et al. (2011), older radial velocity measures are available in the literature, but we decided to consider only the data-sets from Lane et al. (2011), because they are more complete and more radially extended.

In the case of NGC 5139, we merged the largest sample of velocity data available (Reijns et al. 2006) with the sample provided by Mayor et al. (1997), which provides measurements for stars located in the central region of the cluster. A delicate issue regarding this cluster is the controversial position of its center, which plays an important role also in our analysis, because we wish to build a radial-dependent velocity dispersion profile, starting from stellar positions expressed in right ascension and declination (for Reijns et al. 2006). To carry out a proper merging of the two data-sets, we have used the position of the center proposed by Mayor et al. (1997); to calculate radial distances, we have followed the procedure described by van de Ven et al. (2006).

4. Results

4.1. Relaxation classes

The values of the dimensionless parameters and the physical scales of the two families of models determined by the photometric and kinematic fits are presented in Table 3 (the fitting procedure is described in Appendix B). Note that, in general, for the f^(ν) models relatively low values of the concentration parameter Ψ are identified, consistent with the comment given at the end of Appendix A.2. Quantitative information about the best-fit models and the observational profiles, such as the number of the photometric and kinematic points, the values of the relevant reduced chi-squared, and the corresponding residuals, are listed in Table 4. To evaluate the goodness of the fits, in Table 5 we compare the values of the reduced photometric and kinematic chi-squared, denoted by $\stackrel{}{{{\mathit{\chi }}_{\mathrm{p}}^{\mathrm{2}}}^{\mathrm{􏽥}}}$ and $\stackrel{}{{{\mathit{\chi }}_{\mathrm{k}}^{\mathrm{2}}}^{\mathrm{􏽥}}}$ respectively, with the two-sided 90% confidence interval (CI), calculated with respect to the reduced χ²-distribution, characterized by the appropriate number of degrees of freedom (see Appendix B.3 for details).

The surface brightness and the line-of-sight velocity dispersion profiles determined by the fit procedure for the models, together with the observed profiles for the clusters in the first, second, and third relaxation class, are shown in Figs. 1–3, respectively. In the panels, solid lines correspond to the best-fit King models and dotted lines to the best-fit f^(ν) models. The vertical solid line marks the position of the King model projected half-light radius, the dotted one the position of the f^(ν) model projected half-light radius. For the surface brightness profiles, the data from Trager et al. (1995) are indicated with circles, the data from other sources with squares. For the velocity dispersion data, the horizontal bars indicate the length of the radial bin in which the data-points have been calculated; they do not have a role in determining the fit. For each data-point the errors are shown as vertical error bars. Note that, even if we insisted on selecting clusters with a reasonable number of data inside R_e, for about half of the clusters, the kinematic profiles are undersampled in their central region.

In the following part of this subsection we will try to give a general assessment of the quality of the fits in the various cases. As a general rule, for a given cluster the family of models yielding the lowest values of the best-fit photometric and kinematic chi-squared is preferred. In practice, the objective properties of the fits are best obtained by checking directly the values provided in Tables 4 and 5 and by inspection of Figs. 1–3. Table 4 shows that generally, for a given cluster, $\stackrel{}{{{\mathit{\chi }}_{\mathrm{p}}^{\mathrm{2}}}^{\mathrm{􏽥}}}\mathit{>}\stackrel{}{{{\mathit{\chi }}_{\mathrm{k}}^{\mathrm{2}}}^{\mathrm{􏽥}}}$, because the photometric profiles are characterized by a larger number of data-points with reported smaller error-bars. Table 5 shows that NGC 6656 is the only cluster for which the values of $\stackrel{}{{{\mathit{\chi }}_{\mathrm{p}}^{\mathrm{2}}}^{\mathrm{􏽥}}}$ and $\stackrel{}{{{\mathit{\chi }}_{\mathrm{k}}^{\mathrm{2}}}^{\mathrm{􏽥}}}$ for the two models are inside the relevant 90% CI. Within the King modeling, four clusters have $\stackrel{}{{{\mathit{\chi }}_{\mathrm{p}}^{\mathrm{2}}}^{\mathrm{􏽥}}}$ inside the 90% CI and nine have $\stackrel{}{{{\mathit{\chi }}_{\mathrm{k}}^{\mathrm{2}}}^{\mathrm{􏽥}}}$ inside the 90% CI. A visual inspection of the fits (in particular, see the kinematic fits for NGC 362, NGC 3201, NGC 6121, NGC 6656, NGC 6809 and the inner photometric profiles of NGC 3201, NGC 6121, NGC 6218, NGC 6656, and NGC 7078) also suggests that some systematic trends are clearly missed by the two families of models used.

Within the class of relaxed globular clusters, for NGC 362 and NGC 7078 it is evident that the King models cannot reproduce the observed surface brightness profiles, and that the f^(ν) models perform better, especially for describing the outer parts of the cluster. As to the observed velocity dispersion profiles, we see that for NGC 7078 the f^(ν) profile is formally more adequate (at the 99.9% confidence level, which, in the following, we denote by CL), while for NGC 362 the King profile is the closer to the observations (at the 90% CL). We should recall that these two globular clusters are flagged in the Harris catalog as post-core-collapse clusters, although the value of the concentration parameter C is smaller than 2.5. Indeed, the observations indicate that there are some processes that cannot be captured by King models (with particular reference to the shallow cusp in the photometric profile of NGC 7078). By looking at the plots in Fig. 1, NGC 104 and NGC 6341 appear to have both observed surface brightness and velocity dispersion profiles well represented by King models (for NGC 6341 the King $\stackrel{}{{{\mathit{\chi }}_{\mathrm{k}}^{\mathrm{2}}}^{\mathrm{􏽥}}}$ falls within the 90% CI), in spite of the fact that, from Table 5, the fits cannot be considered entirely satisfactory. Curiously, NGC 6121 is equally well described by the two models; the values of the reduced chi-squared are slightly lower for the King model, but in practice the quality of the two fits is similar. We notice that for NGC 104 and NGC 6121 it is particularly evident that the last point in the velocity dispersion profile is significantly higher than expected by both models. A partial explanation of this fact could be the very large extent of the radial interval in which it is calculated. Therefore, it is important to obtain more velocity data-points in the outer regions to clarify this issue.

Fig. 1 Fits by King models and anisotropic f(ν) models to the surface brightness profiles and to the line-of-sight velocity dispersion profiles of relaxed globular clusters. In all panels, solid lines correspond to the King-model fits, dotted lines to f(ν)-model fits; the vertical solid line marks the position of the King model projected half-light radius, Re, the dotted one the position of the f(ν) model projected half-light radius, Re. For the surface brightness profiles, the data from Trager et al. (1995) are indicated with circles, the data from other sources (see Sect. 3.1) with squares. For each data-point, errors are shown as vertical error bars; in the case of the velocity dispersion profile, the horizontal bars indicate the length of the radial bin in which the data-points have been calculated and have no role in the fitting procedure (see Appendix B).

Fig. 2 Fits by King models and anisotropic f(ν) models to the surface brightness profiles and to the line-of-sight velocity dispersion profiles of globular clusters in the intermediate relaxation condition, in the same format as in Fig. 1.

Fig. 3 Fits by King models and anisotropic f(ν) models to the surface brightness profiles and to the line-of-sight velocity dispersion profiles of partially relaxed globular clusters, in the same format as in Fig. 1.

Best-fit profiles for globular clusters in the intermediate relaxation class are shown in Fig. 2. For NGC 288, NGC 6218, and NGC 6809 King models provide a better fit to both the surface brightness and the velocity dispersion profile (at the 90% CL in all cases). For the other globular clusters in this relaxation class the results are less sharp. In fact, for NGC 3201 and NGC 6254 the surface brightness profiles are formally better reproduced by f^(ν) models (for the first, at 95% CL), while the corresponding velocity dispersion profiles are formally better described by King models (at the 95% and 90% CL, respectively). For NGC 6656 both the surface brightness and the velocity dispersion profile are approximately equally well reproduced by the two families of models (at the 90% CL in all cases). For NGC 3201, NGC 6656, and NGC 6809 the velocity dispersion profiles have an irregular shape in the central regions: even if King models formally perform better than f^(ν) models, they are unable to reproduce the observations. We tried to choose a different binning for the data and we found that this irregularity does not depend on the way in which the observed velocity dispersion profile is constructed from the available data-set. Clearly, more data are necessary in order to obtain a more convincing description of these systems.

In Fig. 3, we show the best-fit profiles for partially relaxed globular clusters. For these two clusters we see discordant results: for NGC 2419 f^(ν) models are more adequate for describing the data (at the 99.99% CL and 90% CL for the photometric and kinematic fit, respectively), while for NGC 5139 King models provide a better fit to the observed profiles. Even if formally King models perform better in describing the kinematic profiles of NGC 5139, they do not provide a satisfactory description of the kinematics of the central parts of the cluster (see Fig. 3); in this respect, the f^(ν) models give a better representation of the inner kinematics.

To summarize, we found that, as expected, f^(ν) models tend to perform globally better than King models for the least relaxed globular cluster of our sample, NGC 2419. For NGC 2419, the good performance of f^(ν) models might correspond to the partial relaxation condition of the cluster, consistent with the physical picture that motivates the definition of the f^(ν) models, as outlined in the Introduction. In addition, for three clusters in the second relaxation class (NGC 3201, NGC 6254, and NGC 6656), f^(ν) models are competitive with King models. Furthermore, f^(ν) models can describe well a relatively steep central slope of the velocity dispersion profile, even when the corresponding photometric profile is cored, while King models and other isotropic truncated models (such as Wilson models) are unable to reproduce this kinematical behavior. This fact is evident from the kinematic fits for NGC 2419, NGC 5139 and for possibly one cluster in the intermediate relaxation condition (NGC 6218).

As far as the behavior of the photometric profiles at large radii is concerned, we see that, especially for NGC 104, NGC 362, NGC 6254, NGC 6341, and NGC 7078 King models do not provide a good description of the truncation, as noted in a number of previous studies (in particular, see McLaughlin & van der Marel 2005; Jordi & Grebel 2010; Küpper et al. 2010, in which the outermost parts of the surface brightness profiles are appropriately modeled by using N-body simulations). In some cases, the observed profile falls between the King and the f^(ν) profiles. This suggests that truncated f^(ν) models might behave systematically better than King models for describing these stellar systems. To a large extent, the modification by truncation in phase space of the anisotropic (non-truncated) f^(ν) models is complementary to the generalization of the isotropic (truncated) King models to models characterized by anisotropic pressure, that is, the so-called Michie-King models (in which the truncated Maxwellian is associated with the anisotropic factor of the Eddington models; see Michie 1963; Gunn & Griffin 1979). At this stage, it would be interesting to compare the behavior of the outer photometric and kinematic profiles of the two families of anisotropic truncated models, to evaluate the interplay between truncation and pressure anisotropy in the two different cases. Of course, the simple physical picture offered by King, f^(ν), and Michie-King models still suffers from a number of limitations, as discussed in the Introduction. Such a simple picture is bound to fail in the modeling of clusters in which core collapse has taken place. Therefore, it is not surprising that clusters such as NGC 362 and NGC 7078, even if they belong to the first relaxation class, are not described by isotropic King models as well as expected.

In closing, we wish to reiterate that, in general, the kinematic fits are crucial to assess if a model is actually able to describe a given globular cluster. Unfortunately, the observed velocity dispersion profiles are generally less accurate and less reliable, with respect to the surface brightness profiles; not only the outer parts (radii close to the truncation radius), but also the inner parts (inside the half-light radius) are often not well sampled as would be desired. The present study confirms that in the future it would be desirable to acquire new and better kinematic data.

Fig. 4 Comparison between the values of the structural parameters derived from King and f(ν) best-fit models. The top left panel shows, as crosses, the values of half-mass radii, the top right panel the total masses, the bottom left panel the central mass density and the bottom right panel the mass-to-light ratios. In each panel, the x-axis refers to the values calculated with King models, and the y-axis to those calculated with f(ν) models; the dotted diagonal line denotes the identity relation. The intrinsic quantity rM would be best given in parsecs, but it is reported here in arcseconds for easier comparison with the observed profiles and with the projected half-light radius Re represented in the previous figures.

4.2. King models vs. f ^(ν) models

The values of the relevant structural parameters derived from the best-fit models are presented in Table 6. To compare quantitatively the properties of the best-fit models selected in the two families of models for each globular cluster, we can correlate the values of the derived parameters, such as the half-mass radius r_M, the total mass M, the central mass density ρ₀ and the mass-to-light ratio M/L. These correlations are illustrated in Fig. 4.

For the majority of the globular clusters considered in our sample, the values of the half-mass radius from f^(ν) models are larger than those obtained from King models; only for NGC 362 the opposite is true. For NGC 288, NGC 5139 and NGC 6809 the values calculated from the two models basically coincide. By comparing the values of the total mass and of the central mass density calculated with the two families of models, we see that for the whole sample (with few exceptions) the values calculated with f^(ν) models are larger than those calculated with King models; this fact is not at all surprising, since the f^(ν) models are not truncated. As to the mass-to-light ratios, we see that there is not a tight correlation between the values calculated with King and f^(ν) models, the latter being almost always larger. Similar trends are noted also in the structural properties derived from (isotropic) models characterized by a more spatially extended truncation, such as the Wilson models (see McLaughlin & van der Marel 2005).

In Table 7, we list the values of the core and half-mass relaxation times calculated by using the two best-fit dynamical models for each globular cluster; the clusters are listed in order of increasing King core relaxation times, and the separation in three classes of relaxation here adopted is marked with horizontal lines. Figure 5 illustrates the correlation between the values of these parameters.

Fig. 5 Comparison between the values of the relaxation times derived from King and f(ν) best-fit models. The left panel shows the values of core relaxation times, the right panel the half-mass relaxation times. The format is the same as in Fig. 4.

When considering the core relaxation times (calculated according to Eq. (10) of Djorgovski (1993), as in the Harris 2010, catalog), we see that the original division in the relaxation classes proposed in this paper on the basis of the values listed in the Harris catalog is confirmed. The only exception is NGC 6341: according to the value of the core relaxation time calculated with its best-fit f^(ν) model, this cluster should belong to the second class, rather than to the first. When we list globular clusters according to increasing core relaxation times estimated with the models identified in this paper, we see that the order of relaxed globular clusters changes with respect to that of Table 1. However, there is a general agreement between values of these quantities calculated with the two families of models.

To calculate the half-mass relaxation times we followed the definition of Eq. (5) in Spitzer & Hart (1971), which is based on the half-mass radius. We see that the values from f^(ν) models are larger than those from King models, except for NGC 362 and NGC 6809 (for the latter cluster, the estimated values are approximately equal). In contrast with the case of the core relaxation times, the values of the half-mass relaxation time turn out to be more model-dependent. Such dependence is likely to be due to the differences between the density profiles of the two families of models, which, in general, are less significant in the central regions and become more evident at radii larger than the half-mass radius. The introduction of a truncation for f^(ν) models would also lead to different values of these parameters. We notice that usually the half-mass relaxation time is calculated by inserting in the relevant definition directly the (projected) half-light radius (see, for example, Harris 2010; McLaughlin & van der Marel 2005). By following the same procedure, we found that the values of the half-mass relaxation times are less model dependent, and closer to the values listed in the Harris catalog.

4.3. Comparison with previous studies

In the case of King models, it is interesting to compare the values of the structural parameters found in the present paper with those obtained in previous studies, as summarized in Table 8. By combining the formal errors from both analyses, the comparison of this paper with McLaughlin & van der Marel (2005) is the most significant because their values result from a fitting procedure which is similar to the one we have followed; for this reason, it is also, to some extent, the most surprising. We recall that in this paper the errors on the parameters can be found in Table 3. McLaughlin & van der Marel (2005) list the formal errors on the parameters in their Tables 10 and 12. No similarly detailed comparison with other papers could be made, because in general error analysis is not provided.

We notice that for ten globular clusters (all but NGC 2419, NGC 6254, and NGC 7078) our values of Ψ agree, within the errors, with the values of McLaughlin & van der Marel (2005); for the scale radius, agreement is found for eight objects (all but NGC 2419, NGC 6121, NGC 6254, NGC 6341, and NGC 7078). We may argue that the discrepancy is partly due to the fact that our surface brightness profiles usually contain a larger number of data-points with respect to those of the cited paper (see Sect. 3.1), even when we do not merge the original profiles by Trager et al. (1995) with those of other sources.

The discrepancies between the values of the central surface brightness derived in this paper and those resulting from previous studies are primarily due the correction for the extinction, usually neglected, introduced in our analysis, following McLaughlin & van der Marel (2005). In fact, such discrepancies are not severe; only three objects (NGC 6341, NGC 6809, and NGC 7078) have values which are not consistent, within the errors, with respect to values determined in previous studies. The reason for this is the fact that we added to these profiles the more recent and more accurate data from Noyola & Gebhardt (2006).

The structural parameter for which the differences between the values obtained in our analysis and those in the literature are the most relevant is the central line-of-sight velocity dispersion. For NGC 288, NGC 2419, NGC 3201, NGC 6121, NGC 6254, and NGC 7078 we see that our values agrees, within the errors, with at least one of the values found in the literature (only four of them with McLaughlin & van der Marel 2005). This discrepancy is not at all surprising, because the kinematic data are the most uncertain. This fact gives one further argument for the need for more numerous and more accurate kinematic data for these systems.

4.4. Central slope of the photometric and kinematic profiles

When considering the surface brightness profiles of our selected globular clusters, we note that there are cases in which the observed photometric profiles deviate from the calculated ones at small radii. In particular, by focusing on the innermost regions of the profiles, we see that this is the case for NGC 6121, NGC 6218, and NGC 7078, for which the models are underluminous, and NGC 288, NGC 6656, and NGC 6809, in which the models are overluminous, and NGC 3201, for which the observed central surface brightness appears to oscillate; in spite of these local discrepancies, the global values of the statistical indicators may be satisfactory (see Table 5).

As far as the velocity dispersion profiles are concerned, in four cases (NGC 288, NGC 3201, NGC 6121, and NGC 6656), both models overpredict the central data-points, while, as mentioned in Sect. 4.1, four clusters (NGC 2419, NGC 5139, and, possibly, NGC 6218 and NGC 6254) show a relatively large gradient of the profile in the central regions. At variance with the study of elliptical galaxies, the kinematic profiles of globular clusters are often undersampled inside the half-light radius.

Recently, the central cusps in the observed photometric and kinematic profiles of some globular clusters have been interpreted as clues of the presence of an intermediate-mass black hole (IMBH) in the center of the system (see the analytical model by Bahcall & Wolf 1976, and the N-body simulations by Baumgardt et al. 2005; Noyola & Baumgardt 2011), and a variety of dynamical models (either defined from distribution functions or as solutions of the Jeans equations), have been used in order to constrain the mass of such central object (the best known example is the controversial case of ω Cen, studied by Noyola et al. 2008; van der Marel & Anderson 2010, with different conclusions). However, Vesperini & Trenti (2010) showed that these shallow photometric cusps are not decisive signatures of the presence of an IMBH, and that they might be related to other dynamical processes; moreover, the authors emphasize the fact that the typical accuracy in the data may be insufficient to characterize the slope of the profile as desired.

As shown in the previous sections, the presence of radially-biased pressure anisotropy (which occurs in the context of the family of f^(ν) models or of the Michie-King models, see Michie 1963) can also produce a relatively rapid decline in the central part of the velocity dispersion profile. Therefore, here we are reiterating a point already noted in the literature, that the presence of a central IMBH should not be considered as the only physical explanation of the existence of a central kinematical peak⁵. Unfortunately, our independent conclusion only confirms that the interpretation of this interesting kinematical feature is more model dependent than desired. We recall that f^(ν) models are characterized by a “realistic” anisotropy profile (see Fig. 6 in Trenti & Bertin 2005): the central regions are more isotropic than the outer ones in velocity space, because the models represent a scenario in which violent relaxation has acted more efficiently in the center. This differential kinematical feature is not always present in dynamical models based on the Jeans approach, which, to obtain a fast decline in the velocity dispersion profile in the absence of a central IMBH, usually requires very high values of the anisotropy parameter, even in the central regions of the cluster (see, for example, Sect 4.2 in Noyola et al. 2008; or Sect 5.3 in Lützgendorf et al. 2011).

4.5. Radial orbit instability

Systems in which there is a great amount of radial kinetic energy with respect to tangential kinetic energy are subject to the radial orbit instability.

Polyachenko & Shukhman (1981) (see also Fridman & Polyachenko 1984) introduced a parameter defined as: $\mathit{\kappa }\mathrm{=}\frac{\mathrm{2}{\mathit{K}}_{\mathrm{r}}}{{\mathit{K}}_{\mathrm{T}}}$(1)where K_r and K_T represent the radial and the tangential component of the total kinetic energy, respectively, and argued that when: $\mathit{\kappa }\mathit{>}\mathrm{1.7}\mathrm{\pm }\mathrm{0.25}$(2)radial orbit instability occurs. Actually, different families of models generally have different threshold values for the instability. A detailed discussion of the onset of the radial orbit instability for the family of f^(ν) models can be found in Trenti & Bertin (2005), Trenti et al. (2005), and Trenti & Bertin (2006), where the validity of the criterion expressed in Eq. (2) is discussed.

We calculated the values of the parameter κ for the f^(ν) models that provide the best fit to the observed profiles of the globular clusters of the sample (see Table 9). It is interesting to note that globular clusters belonging to a given relaxation class tend to have similar values of the stability parameter (κ ≈ 1.3,1.8,1.7 for the first, second, and third relaxation class, respectively; in this respect, NGC 6656 and NGC 6254 appear to exhibit an exceptional behavior). In other words, the relaxed class is found to be more isotropic by the f^(ν) diagnostics.

The majority of the globular clusters that have κ ≳ 1.7 are likely to be in a condition of marginal instability. The case of NGC 6254 does remain problematic.

We argue that the introduction of a truncation in phase space to the family of f^(ν) models might have a stabilizing effect, since such truncation will affect primarily the outer parts of a given configuration, which are dominated by radially-biased pressure anisotropy. Therefore the truncation is likely to reduce the global value of the radial component of the total kinetic energy.

This is one more reason to undertake the study of truncated f^(ν) models for application to globular clusters. We plan it for the immediate future.

5. Conclusions

In this paper we have performed a detailed combined photometric and kinematic study of a sample of Galactic globular clusters, representing systems under different relaxation conditions. For these objects, surface brightness and velocity dispersion profiles have been fitted by means of two different families of dynamical models, the truncated isotropic King models and the non-truncated, anisotropic f^(ν) models. The analysis has been carried out by following the same procedure used in the past to study the dynamics of elliptical galaxies. Each globular cluster is then associated with two best-fit models. The main conclusions can be summarized as follows:

• The expected trend, that King models should perform bet-ter for more relaxed globular clusters, has been checkedto exist but it is not as sharp as anticipated. The two clusters(NGC 104 and NGC 6341) forwhich the global fit by King models is most convincing in-deed belong to the class of relaxed objects. King models tendto offer a good representation of the observed photometric pro-files (as is commonly reported), regardless of the relaxationcondition of the system (but a statement of this kind shouldalso be supported by the relevant statistical indicators; seeTable 5). However, the quality of the fits by Kingmodels to the kinematic profiles remains to be proved, even forrelaxed clusters, because of the few data-points and the large errorbars in the observed profiles. Three clusters for which the Kingmodels appear to be inadequate do not actually come as a surprise:NGC 2419 is the least relaxed cluster of the sampleand NGC 362 and NGC 7078 aresuspected to be post-core-collapse clusters.

• The second expected trend, that less relaxed clusters might exhibit the characteristic signature of incomplete violent relaxation, is also partly present but is not as sharp as might have been hoped for. Some cases indeed point to a significant role of radially-biased pressure anisotropy. The least relaxed cluster, NGC 2419, is well described by the f^(ν) models. For the second least relaxed cluster, NGC 5139, the central shallow cusp in the velocity dispersion profile appears to be well captured by the f^(ν) models. A marginal indication in favor of the f^(ν) models also comes from inspection of the inner kinematic profiles of NGC 6218 and NGC 6254, although these two clusters are not among the least relaxed objects. In contrast, King models and other isotropic models (such as the spherical Wilson models) have difficulty in matching significant velocity gradients inside the half-light radius. Therefore, the partial success of the f^(ν) models suggests that for some globular clusters radially-biased pressure anisotropy may be important. This property could be examined further by means of better spatially-resolved kinematic data in the inner regions (i.e., at radii out to approximately the half-light radius). This result is in line with the conclusions of recent papers: in particular, see Ibata et al. (2011) and references therein, based on the application of King-Michie models.

• In some clusters, regardless of the relaxation condition, some qualitative characteristics of the observed profiles are missed by both families of models considered in this paper. It may be that part of these cases would be resolved by a study in terms of truncated f^(ν) models. But it may also be that other ingredients ignored in the paper, such as rotation (solid-body or differential), play a role.

• The paper demonstrates that the values of some structural parameters, such as the total mass and the half-mass radius, can be significantly model-dependent. In view of the results listed in the previous items, this is a clear warning against an indiscriminate use of structural parameters for globular clusters based on only one family of models (the spherical King models).

• In general, the kinematic fits are crucial to assess if a model is actually suited to describe a given globular cluster. The main issue in testing dynamical models on globular clusters is therefore the general lack of good kinematic data: the data are available for a small fraction of the population of Galactic globular clusters and generally made of a small number of data-points, not well distributed in radius. Surprisingly, the kinematic profile is often not well sampled as desired inside the half-mass radius. As discussed in the paper this is a key region for confronting the performance of different dynamical models. In addition, accurate data in the outermost parts, close to the truncation radius,would touch on other important issues, such as the role of tides. Only after the acquisition of good kinematical profiles will it be possible to address properly the issue of dark matter in globular clusters.

Online material

Appendix A: The dynamical models

A.1. Spherical isotropic King models

The King (1966) one-component, spherical, isotropic models are based on the following distribution function: ${\mathit{f}}_{\mathrm{K}}\mathrm{=}\{\begin{array}{c}\\ \mathit{A}\left[\exp \mathrm{(}\mathrm{-}\mathit{aE}\mathrm{)}\mathrm{-}\exp \mathrm{(}\mathrm{-}\mathit{a}{\mathit{E}}_{\mathrm{0}}\mathrm{)}\right]& \mathit{E}\mathrm{\le }{\mathit{E}}_{\mathrm{0}}\\ \mathrm{0}& \mathit{E}\mathit{>}{\mathit{E}}_{\mathrm{0}}\mathit{,}\end{array}$(A.1)

where A, a, E₀ are positive constants, defining two scales and one dimensionless parameter, and E represents the specific energy E = v²/2 + Φ(r), where Φ(r) is the mean-field gravitational potential, to be determined from the Poisson equation. The quantity E₀ is a threshold energy, above which the stars are considered unbound; it can be translated into a truncation radius, r_tr, for the system.

We recall that the radial scale r₀ can be defined as ${\mathit{r}}_{\mathrm{0}}\mathrm{=}{\left(\frac{\mathrm{9}}{\mathrm{4}\mathit{\pi Ga}{\mathit{\rho }}_{\mathrm{0}}}\right)}^{\mathrm{1}\mathit{/}\mathrm{2}}\mathit{,}$(A.2)where ρ₀ is the central mass density. The other dimensional scale is the total mass of the cluster or, as an alternative, its central velocity dispersion.

The dimensionless concentration parameter is given by the central dimensionless potential Ψ = a [Φ(r_tr) − Φ(0)]  or, alternatively, by the index $\mathit{c}\mathrm{=}\log \left(\frac{{\mathit{r}}_{\mathrm{tr}}}{{\mathit{r}}_{\mathrm{0}}}\right)\mathit{,}$(A.3)because these two quantities are in a one-to-one relation.

The quantity listed in the Harris catalog as concentration parameter, in this paper indicated with C, is defined using the standard core radius R_c in place of r₀ in Eq. (A.3). The values of C for a model identified by Ψ ≳ 4 are slightly larger than the corresponding values of c; in fact, for these models 0.8 ≲ R_c/r₀ ≲ 1.

A.2. Anisotropic f^(ν) models

Several families of dynamical models have been developed to represent the final state of numerical simulations of the violent relaxation process thought to be associated with the formation of bright elliptical galaxies via collisionless collapse (for a review, see Bertin & Stiavelli 1993). These models show a characteristic anisotropy profile, with an inner isotropic core and an outer envelope that becomes dominated by radially-biased anisotropic pressure. They provide a good representation of the photometric and kinematic properties of elliptical galaxies. Here we will refer to the family of spherical, anisotropic, non-truncated f^(ν) models (which have been revisited recently in detail by Bertin & Trenti 2003).

The distribution function that defines these models depends on specific energy E and angular momentum J: ${\mathit{f}}^{\mathrm{\left(}\mathit{\nu }\mathrm{\right)}}\mathrm{=}\{\begin{array}{c}\\ \mathit{A}\exp \left[\mathrm{-}\mathit{aE}\mathrm{-}\mathit{d}{\left(\frac{{\mathit{J}}^{\mathrm{2}}}{\mathrm{|}\mathit{E}{\mathrm{|}}^{\mathrm{3}\mathit{/}\mathrm{2}}}\right)}^{\mathit{\nu }\mathit{/}\mathrm{2}}\right]& \mathit{E}\mathrm{\le }\mathrm{0}\\ \mathrm{0}& \mathit{E}\mathit{>}\mathrm{0}\mathit{,}\end{array}$(A.4)

where A, a, d, and ν are positive constants, defining two scales and two dimensionless parameters. For applications, as described by Bertin & Trenti (2003), the dimensionless parameter ν can be fixed at ν = 1. Therefore, similarly to the King models, after integration of the relevant Poisson equation the f^(ν) models are a one-parameter family of models, parametrized by their central concentration, which can be expressed by the central dimensionless potential Ψ =  −aΦ(0).

The physical scales can be expressed as: $\begin{array}{c}\\ {\mathit{r}}_{\mathrm{scale}}\mathrm{=}{\mathit{d}}^{\mathrm{-}\mathrm{1}\mathit{/}\mathit{\nu }}{\mathit{a}}^{\mathrm{-}\mathrm{1}\mathit{/}\mathrm{4}}\\ {\mathit{M}}_{\mathrm{scale}}\mathrm{=}{\mathit{d}}^{\mathrm{-}\mathrm{3}\mathit{/}\mathit{\nu }}{\mathit{a}}^{\mathrm{-}\mathrm{9}\mathit{/}\mathrm{4}}\mathit{.}\end{array}$(A.5)By definition, these models are non-truncated and, because of this, are likely to be less suited to describe the outer parts of globular clusters. A study of globular clusters based on truncated f^(ν) models is postponed to a separate investigation.

We recall that the surface brightness profiles for concentrated models (Ψ ≳ 7) are very close to de Vaucouleurs profile, while for low values of Ψ the models exhibit a sizeable core.

Appendix B: Fitting procedure

In the following we describe the procedure that we have adopted to perform the statistical analysis of the data. We basically follow Bertin et al. (1988).

B.1. Photometric fit

For each family of models considered in this paper, the photometric fit determines which equilibrium model has the projected mass distribution that best reproduces the surface brightness profile, under the assumption that: $\mathrm{\Sigma }\mathrm{\left(}\mathit{R}\mathrm{\right)}\mathrm{=}\frac{\mathit{M}}{\mathit{L}}\mathit{\lambda }\mathrm{\left(}\mathit{R}\mathrm{\right)}\mathit{,}$(B.1)where Σ(R) is the model surface mass density, λ(R) the model surface luminosity density, and the mass-to-light ratio M/L is considered to be constant in the cluster. From the photometric data, we perform a fit that allows us to determine three parameters for each model:

• Ψ, the concentration parameter, which determines the shape of the surface brightness and velocity dispersion profiles;

• r_s, the scale radius, that is r₀ for the King models and r_scale for the f^(ν) models;

• μ₀, the central surface brightness.

We might determine the values of these parameters by minimizing the quantity: $\sum _{\mathit{i}\mathrm{=}\mathrm{1}}^{{\mathit{N}}_{\mathrm{p}}}{\left[\frac{\mathit{m}\mathrm{\left(}{\mathit{R}}_{\mathit{i}}\mathrm{\right)}\mathrm{-}\mathit{S}\mathrm{(}{\mathit{R}}_{\mathit{i}}\mathit{/}{\mathit{r}}_{\mathrm{s}}\mathrm{)}\mathrm{+}{\mathit{\mu }}_{\mathrm{0}}}{\mathit{\delta }{\mathit{m}}_{\mathit{i}}}\right]}^{\mathrm{2}}\mathit{,}$(B.2)where S is the surface density expressed in magnitudes, but we found that it is more convenient to minimize: ${\mathit{\chi }}_{\mathrm{p}}^{\mathrm{2}}\mathrm{=}\sum _{\mathit{i}\mathrm{=}\mathrm{1}}^{{\mathit{N}}_{\mathrm{p}}}{\left[\frac{\mathit{l}\mathrm{\left(}{\mathit{R}}_{\mathit{i}}\mathrm{\right)}\mathrm{-}\mathit{k}\mathrm{\Sigma ̂}\mathrm{(}{\mathit{R}}_{\mathit{i}}\mathit{/}{\mathit{r}}_{\mathrm{s}}\mathrm{)}}{\mathit{\delta }{\mathit{l}}_{\mathit{i}}}\right]}^{\mathrm{2}}\mathit{,}$(B.3)because the chi-squared function turns out to be more stable, as pointed out by McLaughlin & van der Marel (2005). The connection between Eqs. (B.2) and (B.3) is easily made, by taking into account the following relations: $\begin{array}{ccc}& & \mathit{m}\mathrm{\left(}{\mathit{R}}_{\mathit{i}}\mathrm{\right)}\mathrm{-}{\mathit{m}}_{\mathrm{0}}\mathrm{=}\mathrm{-}\mathrm{2.5}\log \mathrm{\left[}\mathit{l}\mathrm{\right(}{\mathit{R}}_{\mathit{i}}\mathrm{\left)}\mathrm{\right]}\\ & & \mathit{S}\mathrm{(}{\mathit{R}}_{\mathit{i}}\mathit{/}{\mathit{r}}_{\mathrm{s}}\mathrm{)}\mathrm{=}\mathrm{-}\mathrm{2.5}\log \mathrm{\left[}\mathrm{\Sigma ̂}\mathrm{\right(}{\mathit{R}}_{\mathit{i}}\mathit{/}{\mathit{r}}_{\mathrm{s}}\mathrm{\left)}\mathrm{\right]}\\ & & \end{array}$(B.4)where the zero-point m₀ = 26.422 allows us to express the observed surface brightness m(R_i) in solar units; the observed luminosities l(R_i) are expressed in solar units and δl_i are calculated from δm_i and represent the errors on luminosities. The quantity is the surface density, normalized to its central value; from this point of view, the parameter k is related to the mass-to-light ratio.

Even if we calculated the parameters from the luminosities, we report the results in terms of magnitudes, because in this way a comparison with the data is more natural. To calculate the surface brightness profile μ(R) that best reproduces the data, we choose the model identified by Ψ, we calculate the normalized projected mass-density profile, and then we rescale it radially with r_s and vertically with μ₀.

Note that we have only two independent fit parameters, because k = k(Ψ,r_s). Hence, the value of the reduced ${\mathit{\chi }}_{\mathrm{p}}^{\mathrm{2}}$ is: $\stackrel{}{{{\mathit{\chi }}_{\mathrm{p}}^{\mathrm{2}}}^{\mathrm{􏽥}}}\mathrm{=}\frac{{\mathit{\chi }}_{\mathrm{p}}^{\mathrm{2}}}{{\mathit{N}}_{\mathrm{p}}\mathrm{-}\mathrm{2}}\mathrm{·}$(B.5)To have a quantitative estimate of the quality of the fit, it is convenient to calculate also the mean and maximum residuals: $\begin{array}{ccc}& & \mathrm{⟨}\mathrm{\Delta }\mathit{\mu }\mathrm{⟩}\mathrm{=}\frac{\mathrm{1}}{{\mathit{N}}_{\mathrm{p}}^{\mathrm{1}\mathit{/}\mathrm{2}}}{\left\{\sum _{\mathit{i}\hspace{0.17em}\mathrm{=}\hspace{0.17em}\mathrm{1}}^{{\mathit{N}}_{\mathrm{p}}}\mathrm{\left[}\mathit{m}\mathrm{\right(}{\mathit{R}}_{\mathit{i}}\mathrm{)}\mathrm{-}\mathit{\mu }\mathrm{(}{\mathit{R}}_{\mathit{i}}\mathrm{)}{\mathrm{]}}^{\mathrm{2}}\right\}}^{\mathrm{1}\mathit{/}\mathrm{2}}\\ & & \end{array}$(B.6)

B.2. Kinematic fit

The parameters identified by the photometric fit determine a model, characterized by the parameter Ψ, the profiles of which are rescaled with r_s. At this point, we perform a fit to the kinematic data to find the central line-of-sight velocity dispersion, that is the velocity scale needed to rescale vertically the normalized projected velocity dispersion profile calculated from the dynamical models. We calculate the value of this scale, V, as the one that minimizes: ${\mathit{\chi }}_{\mathrm{k}}^{\mathrm{2}}\mathrm{=}\sum _{\mathit{i}\hspace{0.17em}\mathrm{=}\hspace{0.17em}\mathrm{1}}^{{\mathit{N}}_{\mathrm{k}}}{\left[\frac{\mathit{\sigma }\mathrm{\left(}{\mathit{R}}_{\mathit{i}}\mathrm{\right)}\mathrm{-}\mathit{V}{\mathit{\sigma }}_{\mathrm{P}}\mathrm{\left(}{\mathit{R}}_{\mathit{i}}\mathrm{\right)}}{\mathit{\delta }{\mathit{\sigma }}_{\mathit{i}}}\right]}^{\mathrm{2}}\mathit{,}$(B.7)where σ_P(R_i) is the projected velocity dispersion calculated from the model, dependent on the rescaled radial coordinate. In this case, the values of the reduced χ² and of the residuals can be found from the following expressions: $\begin{array}{ccc}& & {\begin{array}{}\mathrm{􏽥}\\ {\mathit{\chi }}_{\mathrm{k}}\end{array}}^{\mathrm{2}}\mathrm{=}\frac{{\mathit{\chi }}_{\mathrm{k}}^{\mathrm{2}}}{{\mathit{N}}_{\mathrm{k}}\mathrm{-}\mathrm{1}}\\ & & \mathrm{⟨}\mathrm{\Delta }\mathit{\sigma }\mathrm{⟩}\mathrm{=}\frac{\mathrm{1}}{{\mathit{N}}_{\mathrm{k}}^{\mathrm{1}\mathit{/}\mathrm{2}}}{\left\{\sum _{\mathit{i}\hspace{0.17em}\mathrm{=}\hspace{0.17em}\mathrm{1}}^{{\mathit{N}}_{\mathrm{k}}}\mathrm{\left[}\mathit{\sigma }\mathrm{\right(}{\mathit{R}}_{\mathit{i}}\mathrm{)}\mathrm{-}\mathit{V}{\mathit{\sigma }}_{\mathrm{P}}\mathrm{(}{\mathit{R}}_{\mathit{i}}\mathrm{)}{\mathrm{]}}^{\mathrm{2}}\right\}}^{\mathrm{1}\mathit{/}\mathrm{2}}\\ & & \end{array}$(B.8)One might argue that a more sensible way to perform the fit would be by means of a single fit procedure, using the combined ${\mathit{\chi }}_{\mathrm{tot}}^{\mathrm{2}}\mathrm{=}{\mathit{\chi }}_{\mathrm{p}}^{\mathrm{2}}\mathrm{+}{\mathit{\chi }}_{\mathrm{k}}^{\mathrm{2}}$, defined as the sum of the photometric and of the kinematic contribution. We did perform tests of this different procedure and found that the results are equivalent to those obtained by performing the photometric fit first and then the kinematic fit at fixed Ψ and r_s. This confirms the qualitative expectation that the kinematical data, being less numerous and less accurate with respect to the photometric ones, have little weight in determining Ψ and r_s and are only needed to determine the scale V, that is, the relevant mass-to-light ratio.

Finally, we addressed the issue of whether a different representation of the results of the best-fit models would be more suitable than the one used in Figs. 1–3. Figure B.1 compares the appearance of a fit for a given object using a linear and a logarithmic representation of the radial scale. Based on this and other tests, we preferred to adopt a mixed representation, as used in the main text of the paper.

Fig. B.1 Linear vs. logarithmic representation for the photometric (upper panels) and kinematic (lower panels) profiles of NGC 6121. On the left, profiles are plotted with a linear radial scale, on the right with a logarithmic scale.

Fig. B.2 Confidence regions and confidence intervals on NGC 6121 King model parameters. The top left figure shows the overlap of the confidence regions, corresponding to confidence levels of 68.3% and 95.4% (solid lines), and the curves of constant χ2 (dashed lines); the black dot marks the position of the maximum of the likelihood (that is the minimum of χ2). The top right figure shows the confidence intervals for the likelihood depending on V; the bottom figures show the confidence intervals on the parameter Ψ and r0; the dark grey area corresponds to the confidence level of 68.3%, the light grey area to the confidence level of 95.4%; a solid vertical line marks the position of the maximum of the likelihood (that is the minimum of χ2).

B.3. Goodness of the photometric and kinematic fits

To measure the goodness of the photometric and kinematic fits, we referred to the confidence intervals on the χ²-distribution. The probability density function of the χ² distribution is defined as: $\mathit{f}\mathrm{\left(}\mathit{x}\mathrm{;}\mathit{n}\mathrm{\right)}\mathrm{=}\frac{{\mathrm{e}}^{\mathrm{(}\mathrm{-}\mathit{x}\mathit{/}\mathrm{2}\mathrm{)}}{\mathit{x}}^{\mathrm{(}\mathit{n}\mathit{/}\mathrm{2}\mathrm{-}\mathrm{1}\mathrm{)}}}{{\mathrm{2}}^{\mathit{n}\mathit{/}\mathrm{2}}\mathrm{\Gamma }\mathrm{(}\mathit{n}\mathit{/}\mathrm{2}\mathrm{)}}$(B.9)when x ≥ 0, and zero otherwise. The parameter n corresponds to the number of degrees of freedom and Γ denotes the gamma function (e.g., see Abramowitz & Stegun 1972, Sect. 26.4).

The two-sided confidence interval relative to a confidence level α on the χ² probability density function is defined as the interval $\mathrm{\left[}{\mathit{\chi }}_{\inf }^{\mathrm{2}}\mathit{,}{\mathit{\chi }}_{\sup }^{\mathrm{2}}\mathrm{\right]}$ such that $\mathit{\alpha }\mathrm{=}{\mathrm{\int }}_{{\mathit{\chi }}_{\inf }^{\mathrm{2}}}^{{\mathit{\chi }}_{\sup }^{\mathrm{2}}}\mathit{f}\mathrm{\left(}\mathit{x}\mathrm{;}\mathit{n}\mathrm{\right)}\mathrm{d}\mathit{x}\mathrm{;}$(B.10)and that: $\mathrm{1}\mathrm{-}\frac{\mathit{\alpha }}{\mathrm{2}}\mathrm{=}{\mathrm{\int }}_{\mathrm{0}}^{{\mathit{\chi }}_{\inf }^{\mathrm{2}}}\mathit{f}\mathrm{\left(}\mathit{x}\mathrm{;}\mathit{n}\mathrm{\right)}\hspace{0.17em}\mathrm{d}\mathit{x}\mathrm{=}{\mathrm{\int }}_{{\mathit{\chi }}_{\sup }^{\mathrm{2}}}^{\mathrm{\infty }}\mathit{f}\mathrm{\left(}\mathit{x}\mathrm{;}\mathit{n}\mathrm{\right)}\hspace{0.17em}\mathrm{d}\mathit{x}\mathit{.}$(B.11)

B.4. Errors on the best-fit parameters

For the methods described below, we refer mainly to Press et al. (2007). In the case of the photometric fit, we used the following procedure to calculate the formal errors on the parameters. First we calculated the Hessian matrix: ${\mathit{H}}_{\mathit{ij}}\mathrm{=}\frac{{\mathit{\partial }}^{\mathrm{2}}{\mathit{\chi }}^{\mathrm{2}}}{\mathit{\partial }{\mathit{x}}_{\mathit{i}}\mathit{\partial }{\mathit{x}}_{\mathit{j}}}\mathit{,}$(B.12)where i = 1,2,3 and x₁ = Ψ, x₂ = r_s and x₃ = μ₀; then we calculated the covariance matrix: $\mathrm{\left[}\mathit{E}\mathrm{\right]}\mathrm{=}\mathrm{2}\mathrm{[}\mathit{H}{\mathrm{]}}^{-1}\mathit{.}$(B.13)Finally, we obtained the errors: δx_i = (E_ii)^1/2.

In the case of the kinematic fit, in which ${\mathit{\chi }}_{\mathrm{k}}^{\mathrm{2}}$ depends analytically only on the parameter V, we can immediately calculate: $\mathit{\delta V}\mathrm{=}{\left(\frac{\mathrm{2}}{{\mathit{\partial }}^{\mathrm{2}}{\mathit{\chi }}_{\mathrm{k}}^{\mathrm{2}}\mathit{/}\mathit{\partial }{\mathit{V}}^{\mathrm{2}}}\right)}^{\mathrm{1}\mathit{/}\mathrm{2}}\mathrm{·}$(B.14)In addition to errors, we identified the relevant confidence regions and intervals for the various parameters, after defining a likelihood function Λ = e ^− χ2/2. First, we calculated confidence regions in the (Ψ,r_s) plane, after marginalizing the likelihood function on k, in the following way: $\mathrm{\Lambda }\mathrm{\left(}\mathrm{\Psi }\mathit{,}{\mathit{r}}_{\mathrm{s}}\mathrm{\right)}\mathrm{=}\mathrm{\int }{\mathrm{e}}^{\mathrm{-}{\mathit{\chi }}^{\mathrm{2}}\mathit{/}\mathrm{2}}\mathrm{d}\mathit{k}\mathit{.}$(B.15)We identified the regions corresponding to confidence levels of 68.3% and 95.4%. To check the results, we compared them with the regions defined by the curves of constant χ² with the appropriate values (see Sect. 15.6 in Press et al. 2007). For each globular cluster and for the families of models, we notice that the regions identified with these two methods overlap in a consistent way.

We also calculated the confidence intervals on the parameters Ψ, r_s and V (for the last parameter, the calculation can be done immediately, because the kinematic likelihood depends only on V; for the other parameters, we have to marginalize once more the likelihood function), corresponding to confidence levels of 68.3% and 95.4%. Comparing them with the errors calculated with the covariance matrix method, we found that the results are consistent.

In Fig. B.2 we show an example of the confidence regions and intervals calculated for the King best-fit model for the globular cluster NGC 6121. The top left panel shows the overlap of the confidence regions (solid lines) corresponding to confidence levels of 68.3% and 95.4%, and the curves of constant χ² (dashed lines). It is clear that there is a good agreement between the two relevant pairs of curves. The bottom panels show the confidence intervals calculated on the marginalized likelihood for the parameters Ψ and r₀, and the top right panel those for the parameter V; in these figures the confidence intervals are

shown as shaded areas. By comparing the extent of the dark grey area with the values of the uncertainties on the corresponding parameters (see Table 3), we see that the different methods lead to consistent results.

B.5. Derived parameters

Once the best-fit parameters have been secured, we may calculate other quantities, which represent some important characteristics of globular clusters.

In particular, for each family of models we calculated the core radius R_c, that is the radial position where the surface brightness equals half its central value, the half-mass radius r_M, the total cluster mass M, the central mass density ρ₀, the mass-to-light ratio M/L, the core relaxation time log T_c and the half-mass relaxation time, log T_M. In addition, for King models we calculated the concentration parameter c and the truncation radius r_tr. In turn, for the f^(ν) models we calculated the anisotropy radius r_α, defined as the radial position where α(r_α) = 1, with $\mathit{\alpha }\mathrm{=}\mathrm{2}\mathrm{-}\mathrm{(}\mathrm{⟨}{{\mathit{v}}_{\mathit{\theta }}^{\mathrm{2}}}^{\mathrm{⟩}}\mathrm{+}{\mathrm{⟨}{{\mathit{v}}_{\mathit{\phi }}^{\mathrm{2}}}^{\mathrm{⟩}}}^{\mathrm{)}}/\mathrm{⟨}{{\mathit{v}}_{\mathit{r}}^{\mathrm{2}}}^{\mathrm{⟩}}\mathit{.}$(B.16)The ratio of this radius to the half-mass radius measures how large is the radially-biased anisotropic part of the globular cluster: the smaller r_α, the larger the anisotropic region of the cluster.

Acknowledgments

We are grateful to D.C. Heggie for many constructive remarks that have helped improve the quality of the paper. We would like to thank E. Vesperini for a critical reading of the manuscript and valuable suggestions, M. Bellazzini, P. Bianchini, M. Lombardi, and C. Nipoti for useful comments and conversations, E. Noyola for providing us with the data for ω Cen, and S. E. Arena and M. Trenti for the routines to compute the profiles of the f^(ν) family of models. This work was partially supported by the Italian MIUR.

References

All Tables

All Figures

Fig. 1 Fits by King models and anisotropic f(ν) models to the surface brightness profiles and to the line-of-sight velocity dispersion profiles of relaxed globular clusters. In all panels, solid lines correspond to the King-model fits, dotted lines to f(ν)-model fits; the vertical solid line marks the position of the King model projected half-light radius, Re, the dotted one the position of the f(ν) model projected half-light radius, Re. For the surface brightness profiles, the data from Trager et al. (1995) are indicated with circles, the data from other sources (see Sect. 3.1) with squares. For each data-point, errors are shown as vertical error bars; in the case of the velocity dispersion profile, the horizontal bars indicate the length of the radial bin in which the data-points have been calculated and have no role in the fitting procedure (see Appendix B).

In the text

	Fig. 2 Fits by King models and anisotropic f(ν) models to the surface brightness profiles and to the line-of-sight velocity dispersion profiles of globular clusters in the intermediate relaxation condition, in the same format as in Fig. 1.
In the text

	Fig. 3 Fits by King models and anisotropic f(ν) models to the surface brightness profiles and to the line-of-sight velocity dispersion profiles of partially relaxed globular clusters, in the same format as in Fig. 1.
In the text

Fig. 4 Comparison between the values of the structural parameters derived from King and f(ν) best-fit models. The top left panel shows, as crosses, the values of half-mass radii, the top right panel the total masses, the bottom left panel the central mass density and the bottom right panel the mass-to-light ratios. In each panel, the x-axis refers to the values calculated with King models, and the y-axis to those calculated with f(ν) models; the dotted diagonal line denotes the identity relation. The intrinsic quantity rM would be best given in parsecs, but it is reported here in arcseconds for easier comparison with the observed profiles and with the projected half-light radius Re represented in the previous figures.

In the text

	Fig. 5 Comparison between the values of the relaxation times derived from King and f(ν) best-fit models. The left panel shows the values of core relaxation times, the right panel the half-mass relaxation times. The format is the same as in Fig. 4.
In the text

	Fig. B.1 Linear vs. logarithmic representation for the photometric (upper panels) and kinematic (lower panels) profiles of NGC 6121. On the left, profiles are plotted with a linear radial scale, on the right with a logarithmic scale.
In the text

Fig. B.2 Confidence regions and confidence intervals on NGC 6121 King model parameters. The top left figure shows the overlap of the confidence regions, corresponding to confidence levels of 68.3% and 95.4% (solid lines), and the curves of constant χ2 (dashed lines); the black dot marks the position of the maximum of the likelihood (that is the minimum of χ2). The top right figure shows the confidence intervals for the likelihood depending on V; the bottom figures show the confidence intervals on the parameter Ψ and r0; the dark grey area corresponds to the confidence level of 68.3%, the light grey area to the confidence level of 95.4%; a solid vertical line marks the position of the maximum of the likelihood (that is the minimum of χ2).

In the text