text:
# Round Table Summary: Stellar interferometry as a tool to investigate atmospheres and to compare observations with models

## 1. Introduction

Long-baseline interferometry at optical and near-infrared wavelengths is an emerging technology which is quickly becoming a useful tool to investigate stellar atmospheres and to compare observations with models. A round table on this topic was held on June 20, 2002, during the IAU Symposium 210 ”Modelling of stellar atmospheres”. The goal of this round table was to discuss and stimulate interactions between observations obtained by optical interferometers and theoretical modelling in the field of stellar atmospheres. Here, I summarize aspects that were discussed at the round table, or during the course of the symposium. This summary is based on the individual contributions which can be found in this volume, supplemented by some further literature that was cited at the symposium, without any claim of completeness. This round table summary is, hence, not intended to be a complete review of the field. Optical interferometry, including interferometric techniques as well as scientific studies, have recently been reviewed by Shao & Colavita (1992), Quirrenbach (2001), and Baldwin & Haniff (2002). OLBIN (optical long baseline interferometry news, edited by P. Lawson, and available at olbin.jpl.nasa.gov) is an excellent source of up-to-date information in the field of optical interferometry. Circumstellar envelopes are also studied by maser observations using radio interferometry (Humphreys, this volume). However, radio interferometry was not part of this round table discussion.

In general, information on stellar atmospheres using optical long-baseline interferometry has so far been derived by comparing model predictions to observations and measurements of

* angular sizes of stars and derived absolute radii and effective temperatures
* center-to-limb intensity variations
* asymmetric shapes of stellar surfaces
* horizontal surface inhomogeneities
* circumstellar matter.

These measurements have often been temporally or spectrally resolved, which leads to additional insights into the atmospheres’ structures and dynamics.

In the following, topics of discussion are summarized in the contexts of (1) interferometric facilities and techniques, (2) concepts to study stellar atmospheres, which are applicable to different object classes, or (3) certain classes of objects. There is, however, overlap between these categories.

## 2. Interferometric techniques and facilities

### Interferometric facilities

Results, modelling of results, or predictions using the following interferometric facilities were presented at the symposium.
Name Location   References in this volume     COAST Cambridge, UK   Freytag & Mizuno-Wiedner, Freytag & Finnson, van Belle et al. (b).   GI2T Calern, France   Dominiciano et al., Jankov et al., Stee & Bittar.   IOTA Mt. Hopkins, USA   van Belle et al. (b).   MARK III Mt. Wilson, USA   Quirrenbach & Aufdenberg.   NPOI Flagstaff, USA   Wittkowski.   PTI Mt. Palomar, USA   van Belle et al. (a,b).   KECK M. Kea, USA   van Belle et al. (a,b).   VLTI Paranal, Chile   Dominiciano et al., Jankov et al., Stee & Bittar, Wittkowski, Wittkowski et al.
Further operating interferometric facilities include
CHARA Mt. Wilson, USA   ISI Mt. Wilson, USA   SUSI Narrabri, Australia
For further information on these facilities, see the references mentioned above, Baldwin & Haniff (2002), Quirrenbach (2001), or OLBIN. The above mentioned interferometric facilities are generally available to the stellar atmosphere community through acceptance of observing proposals and/or collaborations, with the exception of the MARK III, which is not operational anymore. Detailed procedures and priorities for scientific topics differ. For instance, ESO currently accepts scientific proposals on all fields from the international community to be carried out on a shared-risk basis with the VLTI and the near-infrared K-band commissioning instrument VINCI. A first call for proposals was issued for Period 70 (October 2002 to March 2003) and a second call for Period 71 (April 2003 to October 2003). Calls for proposals for the scientific instruments MIDI (mid-infrared) and AMBER (near-infrared) can be expected for 2003/2004. In general, due to the limiting magnitudes, interferometers are mostly used for stellar astrophysics, and hence stellar atmospheres. Additional priorities of large facilities using 8–10m class telescopes (Keck, VLTI, LBT) include planet searches and active galactic nuclei.

### Interferometric observables

The primary observables of an interferometer are the amplitude and phase of the complex visibility function, which is the Fourier transform of the object intensity distribution (van Cittert-Zernike Theorem). In presence of atmospheric turbulence, usually only the squared visibility amplitudes and the triple products are accessible by optical and near-infrared interferometers. A triple product is the product of three complex visibilities corresponding to baselines that form a triangle. The phase of the triple product, the closure phase, is free of atmospheric phase noise (Jennison 1958). Thus, the measurement of triple products is a method to obtain phase information of the object intensity distribution, i.e. allows us to detect deviations from point symmetry. Observational results that were discussed at the symposium make use of squared visibility amplitudes (Quirrenbach & Aufdenberg, van Belle et al. (a/b), Wittkowski, Wittkowski et al.). The COAST images of Betelgeuse (Buscher et al. 1990), modelled by Freytag & Finnson and Freytag & Mizuno-Wiedner, as well as the NPOI limb-darkening measurements presented by Wittkowski make in addition use of triple products. The vast majority of interferometric studies use measurements of relatively large (squared) visibility amplitudes before the first minimum of the visibility function. Since this means that less than one resolution element across the stellar disk is obtained, one parameter, i.e. the diameter of the disk, is derived. These diameter measurements can however be done at different position angles in order to detect asymmetric shapes of the stellar disks (van Belle et al. (a,b)), at different wavelengths (Quirrenbach & Aufdenberg, Wittkowski), or at different time periods to detect for instance stellar pulsations. Interferometric observations of more than one resolution element across the stellar disks are difficult to perform, and hence rare. The long baselines needed to obtain this resolution, i.e. to resolve the stellar disk, also produce very low visibility amplitudes. These low visibility amplitudes correspond to vanishing fringe contrasts. Hence, the corresponding fringes are difficult to be detected, tracked, and recorded. In addition, the data reduction of low visibility values is often a special challenge as well, since the contributions of photon noise and detection bias terms might become dominant. Such observations and the successful compensation of bias terms are discussed in this volume by Wittkowski.

In addition to visibility amplitudes and triple products, the phase of the complex visibility is an interferometric observable as well, as mentioned above. If atmospheric and instrumental noise can be controlled or compensated, for example by phase referencing using a reference star within the isoplanatic patch or by using the same source at another wavelength, the relative position of the object on the sky can be derived with an accuracy which is much higher than the diffraction limit. If the position of the photocenter of the object depends on the wavelength, such differential phase measurements can be used to derive information on the spatial object brightness distribution. Information on stellar diameters, stellar rotation, and starspots can be derived. This technique was described by Beckers (1982), Chelli & Petrov (1995), Jankov et al. (2001), and is discussed in this volume by Dominiciano de Souza et al. and Jankov et al. Jankov et al. discuss in particular the imaging potential of a tomographic technique which combines time-resolved spectroscopy (Doppler Imaging) and differential interferometry.

Simulations of complete atmospheres of stars, such as the radiation hydrodynamics simulations of Betelgeuse presented by Freytag & Finnson and Freytag & Mizuno Wiedner would ideally be compared with high-resolution images of stellar surfaces derived by interferometric observations. If modulus and phase of the complex visibility values are derived, the object intensity distribution can in principle be directly reconstructed from interferometric data using imaging techniques which effectively interpolate the limited coverage of the Fourier-plane. A recent review describing the image fidelity using optical interferometers can be found in Baldwin & Haniff (2002). In order to reconstruct an image of reasonable fidelity, it is intuitively understandable that the number of visibility data points has to be at least as large as the number of unknowns, i.e. the pixel intensities of the image. Furthermore, despite the ability of imaging algorithms to effectively interpolate the sparse aperture data, the aperture has to be filled as uniformly as possible in order to avoid strong artifacts of the PSF such as sidelobes. This means that many very low visibility values at several lobes of the visibility function and at different azimuth angles have to be measured. As a result, images will probably be limited to a few pixels (about 3 $`\times `$ 3 – 5 $`\times `$ 5 px.) in the next years. It was discussed that model fitting, for instance fitting of parameters that describe asymmetries as a function of time, might be a better choice for first attempts.

### Calibration

Interferometric observations are generally limited by the calibration of the atmospheric and instrumental transfer function, i.e. the instrumental visibility response to a point source. This quantity is, dependent on instrument and site, observed to be dependent on different parameters as for instance seeing, atmospheric coherence time, or zenith angle of the source. The transfer function is usually calibrated by observations of an unresolved star or a star with well-known diameter. Since, especially if long baselines are used, almost all calibrator stars are at least slightly resolved, the uncertainty of the knowledge of the diameter of the calibrator star becomes often a major limitation (e.g. Percheron et al. 2002). In addition, the diameter of the calibrator star is often obtained at a wavelength different from that of the new experiment. In this case, model atmospheres need to be used to predict the interferometric transfer function. This can best be done by using the Rosseland mean angular diameter. Then, the model atmosphere is needed to derive the transfer function based on the Rosseland mean angular diameter, or to transform the Rosseland mean angular diameter into a more practical uniform disk diameter at the wavelength of the experiment. As a result, calibration is a topic for which optical interferometry relies on stellar atmosphere modelling, rather than providing observational constraints for models.

## 3. Concepts to study stellar atmospheres

The following concepts to study stellar atmospheres were discussed. These concepts are usually applicable to different classes of objects.

### Diameters, radii, effective temperatures

The vast majority of interferometric measurements have so far directly derived only the first stellar surface structure parameter, which is the angular size of the star. The combination of the angular diameter with the bolometric flux is the most direct method to obtain the effective temperature of the star. This is an important constraint for comparisons of model atmospheres with stellar spectra, since it fixes one free parameter. Another diameter that can be derived from the angular diameter is the absolute radius of the star if its distance is known, e.g. by using the Hipparcos parallax. Further fundamental stellar parameters have been derived from observations of binary systems (e.g. Hummel et al. 1995). The stellar diameter is usually obtained using a uniform disk model. The uniform disk diameter is wavelength-dependent due to the limb-darkening effect (see, e.g. Quirrenbach et al. 1996, Wittkowski et al. 2001). Stellar model atmospheres are used to transform wavelength-dependent uniform disk diameters at the observational wavelength into wavelength-independent limb-darkened diameters. The observed relative variation of uniform disk diameters as a function of wavelength can also be compared to model predictions of this wavelength dependence, and hence be used to test model atmospheres. Quirrenbach & Aufdenberg present in this volume uniform disk diameter measurements of 47 cool giants with the Mark III Interferometer on Mt. Wilson, California in a strong TiO band and in the continuum. They reproduce the diameter ratios by PHOENIX state of the art model atmospheres (Hauschildt et al. 1997; Baron et al., this volume) and derive important constraints for these models. This concept can be used for many different stars using modern interferometers equipped with instruments that allow spectral resolution, as for instance the VLTI with the AMBER instrument at near-infrared wavelengths. Van Belle et al. (b) discuss in this volume wavelength-dependent diameter measurements of Mira stars. Narrowband diameters of oxygen-rich and carbon-rich Miras and non-Mira stars obtained with the PTI are presented; Mira diameter measurements obtained by different interferometers are reviewed.

### Asymmetric shapes of stellar surfaces

Interferometric measurements can also derive asymmetric shapes of stellar surfaces if baselines of different orientation are used. Van Belle et al. (a) present in this volume the recent detection of an oblate photosphere of a fast rotator, the main sequence star Altair (see also van Belle et al. 2001). The theoretical description and simulations of observations of stellar rotation are also discussed in this volume by Dominiciano de Souza et al. and Jankov et al. (see also Domiciano de Souza et al. 2002, Jankov et al. 2001). Van Belle et al. (b) discuss in this volume that also for Mira stars indications of departures from spherical symmetry are beginning to be observed as increasingly rich image information is obtained by a new generation of interferometers. Mira stars rotate slowly, and hence, stellar rotation is an unlikely cause for these asymmetries. The nature of the asymmetry of surfaces of AGB stars is not well understood.

### Intensity profiles and limb-darkening

Optical interferometry has proven its capability to go beyond the measurement of diameters, and to measure additional surface structure parameters. Through the direct measurement of stellar limb-darkening, interferometry tests the wavelength-dependent intensity profile across the stellar disk. Unfortunately, direct measurements of the limb-darkening effect are rare because of the observational difficulties mentioned above. While diameters have so far been obtained for several hundred stars with interferometric and lunar occultation techniques, limb-darkening has been directly observed for only a very limited number of stars (Hanbury Brown et al. 1974, Di Benedetto & Foy 1986, Haniff et al. 1995, Quirrenbach et al. 1996, Burns et al. 1997, Hajian et al. 1998, Wittkowski et al. 2001). Wittkowski (this volume) present multi-wavelength limb-darkening measurements on the giant star $`\gamma `$ Sge, obtained with the NPOI, which succeeded not only in discriminating between uniform disks and limb-darkened disks, but also in constraining Kurucz model atmosphere parameters (Wittkowski et al. 2001). First VLTI measurements of visibility values beyond the first minimum were shown as well. Confronting models with both measured stellar spectra and measured intensity profiles at the same time is a strong test for the radiation fields predicted by atmosphere models (see, e.g. Aufdenberg & Hauschildt 2003). Considering the rapid evolution of interferometry, it can be expected that more limb-darkening parameters of a wider variety of stars can be measured soon with high spectral resolution.

### Surface features

The final goal is the reconstruction of images of stellar surfaces with several resolution elements across the stellar disk, revealing both the overall center-to-limb intensity variation and additional horizontal inhomogeneities as warm or cool starspots. These inhomogeneities may be caused by photospheric convection, magnetic fields, or abundance inhomogneities. Small scale features on stellar surfaces have already been suggested by interferometric observations of Capella (Hummel et al. 1994) and $`\beta `$ And (Di Benedetto & Bonneau 1990). For the apparently largest supergiants $`\alpha `$ Orionis, $`\alpha `$ Scorpii and $`\alpha `$ Herculi, a few bright spots have been detected and mapped by direct imaging techniques including optical interferometry and HST imaging (Buscher et al. 1990, Wilson et al. 1992, Gilliland & Dupree 1996, Burns et al. 1997, Tuthill et al. 1997, Young et al. 2000). It is expected that starspots on relatively large stars can directly be detected and its parameters be constrained using squared visibility amplitudes and triple products obtained with the VLTI with the AMBER instrument (see Wittkowski, this volume; Wittkowski et al. 2002). Jankov et al. discuss in this volume the potential to combine the techniques of stellar differential (phase) interferometry and Doppler imaging in order to constrain surface features on smaller rotating stars. Parameters describing horizontal surface inhomogeneities can be compared to simulations of entire stellar atmospheres. Freytag & Mizuno-Wiedner and Freytag & Finnson (this volume) present 3D radiation hydrodynamic simulations of the envelope and atmosphere of a red supergiant, in order to model observational data of Betelgeuse.

## 4. Classes of objects

A few particular classes of objects were discussed in the course of the symposium, for which comparisons of models with interferometric observational data seem to be promising.

### Red giants and supergiants

Red giant and supergiant stars have so far been the prime target for interferometric observations. This is caused by their large apparent sizes, and their brightness, especially at near-infrared wavelengths. Surface features have been detected on the surfaces of the apparently largest supergiants $`\alpha `$ Orionis, $`\alpha `$ Scorpii and $`\alpha `$ Herculi (e.g. Buscher et al. 1990, Burns et al. 1997, Young et al. 2000). Models to describe such surface features are presented by Freytag & Mizuno-Wiedner and Freytag & Finnson (this volume).

### AGB stars and Mira stars

The evolution of late-type stars along the asymptotic giant branch (AGB) is characterized by high luminosities and low effective temperatures, and accompanied by significant mass loss to the circumstellar environment with mass-loss rates of up to $`10^710^4`$ M/year. AGB stars with masses of the order of 1 M become unstable to large amplitude radial pulsations and become Mira variable stars. The processes of mass-loss and pulsation are essential for our understanding of late stages of stellar evolution, but are still a matter of debate (see the excellent keynote presentation by Woitke in this volume). Atmosphere models need to couple time-dependent dynamics, radiative transfer, chemistry, and dust formation (Woitke, Höfner et al., both in this volume). Optical interferometry can provide important observational constraints by measuring the intensity profile across the stellar disks, which are often observed to be very different from uniform disks (e.g., Wittkowski et al., this volume), the time/stellar phase dependent diameter variations, and absolute radii indicating the mode of pulsation (see, e.g. Hofmann et al. 2002). Interferometric observations of Mira stars have recently been performed by COAST (Young et al. 2000), IOTA (van Belle et al. 2002, Hofmann et al. 2002), and PTI (Thompson et al. 2002). These observations are reviewed by van Belle et al. in this volume. Wittkowski et al. present first preliminary results of Mira observations with the VLTI and the VINCI commissioning instrument.

### Small cool stars

Ludwig (this volume) present model atmospheres for mid M-type main-, as well as pre-main-sequence objects. These model atmospheres are characterized by low effective temperatures and high surface gravities. Because of the rather small angular sizes, these objects are difficult for interferometric observations. However, first measurements of diameters of M dwarfs succeeded already (Lane et al. 2001, Segransan et al. 2003). These observations provide an empirical mass-radius relation, and can thus provide observational constraints for atmosphere models of this parameter space (low masses, low effective temperatures, high surface gravities).

### Winds from Be stars

Classical Be stars are rapidly rotating hot giant or main-sequence stars. They show asymmetric circumstellar winds in the light of emission lines. Stee & Bittar (2001, and this volume) present theoretical models of Be star winds, and predict constraints by near-infrared and visible interferometric observations. Optical interferometry has already been successfully used to study Be star winds (Thom et al. 1986, Quirrenbach et al. 1993, Stee et al. 1995). Upcoming facilities with high spatial and high spectral resolution as for instance VLTI with the AMBER instrument are expected to provide stronger constraints on Be star winds.

## 5. Concluding remarks

Stellar atmosphere models are still mostly constrained by observations of stellar spectra only, which are integrated over the stellar disks. This round table discussion has shown that optical interferometry has already proven its ability to provide important additional observational constraints by spatially resolved observations of stellar disks and circumstellar matter. The comparison of model atmospheres with both stellar spectra and interferometric observations is a strong test and leads as well to increased confidence in these models. Observational results for different object classes and for different parameter spaces of stellar model atmospheres have been presented in the course of this symposium. Optical interferometry is currently rapidly evolving, and simulations show that stronger constraints for more stars can be expected soon. Most progress in this field can be made by increasing interactions between the efforts of theoretical atmosphere modelling and interferometric observations.

### Acknowledgments.

I am grateful to the speakers of this round table discussion, Bernd Freytag, Robert Kurucz, Andreas Quirrenbach, Slobodan Jankov, Philippe Stee, Gerard van Belle, for their contributions. I thank the SOC for allocating time for a round table on stellar interferometry, and for the invitation as a round table organizer.

## References

Aufdenberg J., Hauschildt P.H., 2003, Proc. SPIE 4838, in press

Baldwin J.E., Haniff C.A., 2002, Phil. Trans. A, 360, 969

Baron E., Hauschildt P.H., Allard F., et al., this volume

Beckers J.M., 1982, Optica Acta 29, 361

Burns D., Baldwin J.E., Boysen R.C., et al., 1997, MNRAS 245, 7

Buscher D.F. Haniff C.A., Baldwin J.E., Warner P.J., 1990, MNRAS 245, 7p

Di Benedetto G.P., Foy R., 1986, A&A 166, 204

Di Benedetto G.P., Bonneau D., 1990, ApJ 358, 617

Dominiciano de Souza A., Vakili F., Jankov S., et al., 2002, A&A 393, 345

Dominiciano de Souza A., Zorec J., Jankov S., Vakili F., Petrov R.G., Abe L., this volume

Freytag B., Finnson S., this volume

Freytag B., Mizuno-Wiedner M., this volume

Gilliland R.L., Dupree A.K., 1996, ApJ 463, L29

Hajian A.R., Armstrong J.T., Hummel C.A., et al., 1998, ApJ 496, 484

Hanbury Brown R., Davis J., Lake R.J.W., et al. 1974, MNRAS 167, 475

Haniff C.A., Scholz M., Tuthill P.G., 1995, MNRAS 276, 640

Hauschildt P.H., Baron E., Allard F., 1997, ApJ 483, 390

Hofmann K.-H., et al., 2002, New Astronomy 7, 9

Höfner S., Sandin C., Aringer B., et al., this volume

Hummel C.A., Armstrong J.T., Quirrenbach A., et al., 1994, AJ 107, 1859

Hummel C.A., Armstrong J.T., Buscher D.F., et al., 1995, AJ 110, 376

Humphreys E.M.L., 2001, this volume

Jankov S., Vakili F., Domiciano de Souza A., this volume

Jankov S., Vakili F. Dominiciano de Souza Jr., Janot-Pacheco E., 2001, A&A 377, 721

Jennison R.C., 1958, MNRAS 118, 276

Kurucz R., 1993, Kurucz CD ROMs, Cambridge, Mass.: Smithsonian Astrophysical Observatory

Lane B.F., Boden A.F., Kulkarni S.R., 2001, ApJ 551, L81

Ludwig H.-G., this volume

Percheron I., Richichi A., Wittkowski M., 2003, Proc. SPIE 4838, in press

Quirrenbach A., Mozurkewich D., Armstrong T., et al., 1992, A&A 259, L19

Quirrenbach A. Mozurkewich D., Armstrong T., et al., 1993, ApJ 406, 215

Quirrenbach A., Hummel C.A., Buscher D.F., et al., 1993, ApJ 416, L25

Quirrenbach A., Mozurkewich D., Buscher D.F., et al., 1996, A&A 312, 160

Quirrenbach A., 2001, Ann. Rev. Astron. Astrophys. 39, 353

Schwarzschild M., 1975, ApJ 195, 137

Segransan D., Kervella P., Forveille T., Queloz D., 2003, A&A, in press

Shao M., Colavita M.M., 1992, Ann. Rev. Astron. Astrophys. 30, 457

Stee P., de Araujo F.X., Vakili F., et al., 1995, A&A 300, 219

Stee Ph., Bittar J., 2001, A&A 367, 532

Thompson R.R., Creech-Eakman M.J., van Belle G.T., 2002, ApJ 577, 447

Tuthill P.G., Haniff C.A., Baldwin J.E., 1997, MNRAS 285, 529

van Belle G.T., Ciardi D.R., Thompson R.R., et al., 2001, ApJ 559, 1155

van Belle G.T., Thompson R.R., Creech-Eakman M.J., 2002, AJ 124, 1706

van Belle G.T., Ciardi D.R., Thompson R.R., Akeson R.L., this volume (a)

van Belle G.T., Thompson R.R., Creech-Eakman M.J., this volume (b)

Wilson, R.W., Baldwin, J.E., Buscher, D.F., Warner, P.J., 1992, MNRAS 257, 369

Wittkowski M, Hummel C.A., Johnston K.J., et al., 2001, A&A 377, 981

Wittkowski M., Schöller M., Hubrig S., et al., 2002, AN 323, 241

Wittkowski M., this volume

Wittkowski M., Richichi A., Correia S., Schöller M., this volume

Woitke P., this volume

Young J.S., Baldwin J.E., Boysen R.C., et al., 2000, MNRAS 315, 625

Young J.S., Baldwin J.E., Boysen, R.C., et al., 2000, MNRAS 318, 381


format: text

text:
# On the linear term in the nuclear symmetry energy

## Abstract

The nuclear symmetry energy calculated in the RPA from the pairing plus symmetry force Hamiltonian with equidistant single-nucleon levels is for mass number $`A=48`$ approximately proportional to $`T(T+1.025)`$, where $`T`$ is the isospin quantum number. An isovector character of the pair field assumed to produce the observed odd-even mass staggering is essential for this result. The RPA contribution to the symmetry energy cannot be simulated by adding to the Hartree-Fock-Bogolyubov energy a term proportional to the isospin variance in the Bogolyubov quasiparticle vacuum, and there are significant corrections to the approximation which consist in adding half the isocranking angular velocity. The present calculation employs a smaller single-nucleon level spacing than used in a previous investigation of the model.

The so-called symmetry energy in the semi-empirical formula for nuclear masses is quadratic in the difference $`NZ`$ between the numbers of neutrons and protons . In a more general view including isobaric analogue states, this may be seen as a quadratic dependence on the isospin quantum number $`T`$. Nuclei with small values of $`T`$ have a larger binding energy than expected from a quadratic extrapolation from higher $`T`$-values . In one parametrization of this extra binding energy, the so-called Wigner energy, a term linear in $`T`$ is included in the total energy so that the symmetry energy becomes proportional to $`T(T+x)`$ for some constant $`x`$. The best choice of $`x`$ seems to be close to 1 , possibly with a preference for a slightly larger value . Other parametrizations of the Wigner energy are employed in Refs.  and .

The origin of the Wigner energy is a debated issue . For a schematic, isobarically invariant Hamiltonian, I recently showed that the Hartree-Bogolyubov energy rises quadratically in $`T`$, whereas the $`T`$-dependence of the RPA correlation energy is dominated by a term equal to half the derivative of the Hartree-Bogolyubov energy with respect to $`T`$. The Hartree-Bogolyubov energy and the dominant $`T`$-dependent term in the RPA correlation energy thus give a symmetry energy proportional to $`T(T+1)`$ .

For the details of the theory, see Ref. . The Hamiltonian is the sum of an independent-nucleon Hamiltonian with equidistant, fourfold degenerate single-nucleon levels, an isobarically invariant isovector pairing force and a symmetry force:

$$H=H_0G𝐏^{}𝐏+\frac{\kappa 𝐓^2}{2},H_0=\underset{k\sigma \tau }{}ϵ_ka_{k\sigma \tau }^{}a_{k\sigma \tau }^{}.$$

Here, $`𝐏`$ is a pair annihilation isovector, $`𝐓`$ is the isospin, and $`G`$ and $`\kappa `$ are coupling constants. The single-nucleon energy $`ϵ_k`$ takes $`\mathrm{\Omega }`$ equidistant values separated by $`\eta `$, and $`a_{k\sigma \tau }`$ are single-nucleon annihilation operators. $`k\sigma `$ is $`k`$ or $`\overline{k}`$, where the bar denotes time reversal, and $`\tau `$ distinguishes neutrons from protons. Hartree-Bogolyubov optimal states are derived from a Routhian

$$R=H\lambda A_\text{v}\mu T_z,$$

where $`A_\text{v}`$ is the number of valence nucleons. The nucleon chemical potential $`\lambda `$ is placed midway between the lowest and the highest $`ϵ_k`$, so that $`A_\text{v}=2\mathrm{\Omega }`$ in the optimal state, and the isocranking angular velocity $`\mu `$ varied to produce a range of $`T_z`$.

The symmetry of even $`A_\text{v}/2\pm T_z`$ is imposed on the trial Bogolyubov quasiparticle vacuum. The model thus describes isobaric multiplets whose substate with $`T_z=T`$ belongs to a doubly even nucleus. Due to this symmetry, the $`z`$-coordinate of the gap isovector $`𝚫=G𝐏`$ vanishes. By a suitable choice of phases, the optimal state is then an ordinary product of neutron and proton BCS states with a common real neutron and proton gap $`\mathrm{\Delta }`$. For numeric convenience, $`\mathrm{\Delta }`$ is kept constant with the variation of $`\mu `$, which results in a negligible variation of $`G`$.

The RPA calculation employs the standard technique of a perturbative boson expansion of products of nucleon field operators . By truncating the resulting expansions of $`H`$ and $`R`$ to second order, RPA operators $`\stackrel{~}{H}`$ and $`\stackrel{~}{R}`$ are obtained. The initial symmetries of $`H`$ and $`R`$ are recovered by replacing in $`\stackrel{~}{H}`$ and $`\stackrel{~}{R}`$ truncated expansions of $`A_\text{v}`$ and $`𝐓`$ by formal variables which obey the mutual commutation relations of the exact operators and have the same commutators with other boson operators as the truncated ones. When $`A_\text{v}`$ and $`T_z`$ are equal to eigenvalues of $`A_\text{v}`$ and $`T_z`$ respectively, the lowest eigenstate of the approximate Routhian thus obtained has $`A_\text{v}=A_\text{v}`$ and $`T=T_z=T_z`$. This state is also an eigenstate of the approximate Hamiltonian with an eigenvalue $`E`$ that is the sum of two terms: the Hartree-Bogolyubov energy

$$E_0=H_0\frac{2\mathrm{\Delta }^2}{G}+\frac{\kappa T^2}{2}$$

and the RPA correlation energy

$$E_2=G\mathrm{\Delta }𝐏^{}\mathrm{\Delta }𝐏+\frac{\kappa \mathrm{\Delta }𝐓^2}{2}+\frac{1}{2}\left(\underset{\nu }{}\omega _\nu \underset{i<j}{}[b_{ij}^{},[\stackrel{~}{R},b_{ij}^{}]]\right).$$
(1)
Here, $`\mathrm{\Delta }𝐏=𝐏𝐏`$ and $`\mathrm{\Delta }𝐓=𝐓𝐓`$, $`\omega _\nu `$ denotes the eigenfrequencies of the boson system with Hamiltonian $`\stackrel{~}{R}`$, and $`b_{ij}^{}`$ are the boson annihilation operators associated with pairs of Bogolyubov quasiparticles. $`E_0`$, $`E_2`$ and $`E=E_0+E_2`$ are defined for any $`T=T_z`$ by these expressions.

In the calculation of Ref. , which was intended to simulate the isobaric chain with mass number $`A=48`$, I used the Fermi gas estimate of the single-nucleon level spacing $`\eta `$. As pointed out by Satula and Wyss  and Głowacz, Satula and Wyss , this gives an unrealistically large value. I have therefore repeated the calculation for two smaller values: $`\eta =4/(6A/(c\text{MeV})/\pi ^2)`$ with $`c=8`$ and 10 , which gives $`\eta =1.1`$ and 1.4 MeV respectively. Otherwise, the parameters of the present calculation are the same as in Ref. : The number $`\mathrm{\Omega }`$ of equidistant, fourfold degenerate single-nucleon levels is equal to 24, the pair gap $`\mathrm{\Delta }=12\text{MeV}/\sqrt{A}=1.7\text{MeV}`$ , and the total coefficient of $`T^2/2`$ in the Hartree-Bogolyubov part of the symmetry energy $`\eta +\kappa =2(134.4A^1203.6A^{4/3})\text{MeV}=3.3\text{MeV}`$ .

The result for the smallest single-nucleon level spacing $`\eta =1.1\text{MeV}`$ is shown in Fig. 1.

A comparison with Fig. 1 of Ref.  shows that the change of $`\eta `$ does not alter the situation qualitatively.111Due to a computation error, $`\mu /2(E_2E_{2,T=0})`$ as shown in Fig. 1 of Ref.  is too large by an almost constant factor about 1.2. Still, almost exactly $`E_0E_{0,T=0}=(\eta +\kappa )T^2/2`$. Since the pairing energy $`2\mathrm{\Delta }^2/G`$ is essentially constant—it varies by less than .1 MeV in the range of $`T`$ considered—this means that even in the presence of pairing, the ‘kinetic’ contribution to the symmetry energy, that is, the contribution from $`H_0`$, is in a very good approximation equal to $`\eta T^2/2`$. This expression is easily shown to be exact for $`\mathrm{\Delta }=0`$.

Also like in the previous calculation, the RPA contribution to the symmetry energy $`E_2E_{2,T=0}`$ is dominated by the frequency $`\mu `$ of the normal mode associated with the direction of the isospin. This is indeed more accurately true here than for larger values of $`\eta `$. Thus, from $`T=0`$ to 4 the sum of all the other terms in $`E_2`$ only change by about .5 MeV for $`\eta =1.1\text{MeV}`$ and about .8 MeV for $`\eta =1.4\text{MeV}`$ as compared to about 1.5 MeV for the value $`\eta =2.1\text{MeV}`$ given by the Fermi gas estimate. The reason for this seems obvious: With the smaller level spacing $`\eta `$, the isodeformation $`\mathrm{\Delta }/\eta `$ is larger and the isorotation therefore more collective. This weakens the coupling between the collective and intrinsic degrees of freedom. It may be noticed that since $`\mathrm{\Delta }A^{1/2}`$ and $`\eta A^1`$, the isodeformation increases with $`A`$. Heavier nuclei should therefore show an even more collective isorotation.

The identification of $`T`$ with $`T_z`$ provides an extension of $`E_2(T)`$ to negative $`T`$. Due to the isobaric invariance of $`H`$, this results in an even function. The only term in $`E_2`$ which is not analytic at $`T=0`$ is $`\mu /2`$, which is replaced by $`\mu /2`$ for negative $`T`$. It follows that $`E_2E_{2,T=0}\mu /2`$ vanishes quadratically for $`T0^+`$. Now, apart from a negligible correction due to the fact that $`\mathrm{\Delta }`$ rather than $`G`$ is kept constant in the calculation, we have $`\mu =dE_0/dT=(\eta +\kappa )T`$. The curve in Fig. 1 is in fact well approximated by the second order polynomial $`E_2E_{2,T=0}=(\eta +\kappa )T/2.033\text{MeV}\times T^2`$. Altogether, we thus have in a good approximation

$$EE_{T=0}=\left(\frac{\eta +\kappa }{2}.033\text{MeV}\right)T(T+x),$$

where

$$x=\left(1\frac{.033\text{MeV}}{(\eta +\kappa )/2}\right)^1=1.020.$$

For $`\eta =1.4\text{MeV}`$, this analysis gives $`x=1.031`$. These values are consistent with the possible indications in the data of $`x`$ being slightly larger than 1 .

From my results in Ref. , Głowacz, Satula and Wyss infer $`x.8`$ for $`T=2`$ . These authors may have calculated $`(E_2E_{2,T=0})T/(E_0E_{0,T=0})`$, which is evidently not a right way of extracting $`x`$ from the energies.

Nuclear masses are extensively compared with Hartree-Fock-Bogolyubov calculations. See for example Ref.  and references therein. If we neglect the change of optimal state by the exchange terms in the self-consistent single-nucleon potential and pair field, the difference between the Hartree-Fock-Bogolyubov and Hartree-Bogolyubov energies, the Fock term, is in the present theory precisely the first term in Expr. (1). This is a somewhat artificial term since it is cancelled by a part of the second sum in the second term.222In Ref. , in the third paragraph of the second colomn on page 289, ‘second term in the second sum’ should be as here ‘second sum in the second term’. It is understood there that if $`\mathrm{\Delta }=0`$, also $`G=0`$. In fact Eq. (1) may be recast into the form

$$E_2=\frac{1}{2}\left(\underset{\nu }{}\omega _\nu \underset{i<j}{}(E_i+E_j)\right),$$
(2)
where $`E_i`$ is the Bogolyubov quasiparticle energy.

Satula and Wyss seem to suggest that the RPA correlation energy may be simulated by adding $`\eta \mathrm{\Delta }𝐓^2/2`$ to the Hartree-Fock-Bogolyubov energy . Neglecting again the change of optimal state, this amounts in the present theory to replacing $`E_2`$ by

$$E_2^{}=G\mathrm{\Delta }𝐏^{}\mathrm{\Delta }𝐏+\frac{(\eta +\kappa )\mathrm{\Delta }𝐓^2}{2}.$$

Since this is a modification of the Fock term, which is absent in Eq. (2), a relationship between $`E_2`$ and $`E_2^{}`$ is not obviously anticipated. Fig. 2

shows that in fact they depend on $`T`$ very differently. In particular, $`E_2^{}(T)`$ is even and analytic at $`T=0`$, so $`E_2^{}E_{2,T=0}^{}`$ vanishes quadratically for $`T0`$. Thus $`E_2^{}`$ does not give the empirical cusp at $`T=0`$ in the graph of the function $`E(T)`$.

The $`T`$-dependence of $`E_2^{}`$ mainly reflects that of the isospin variance $`\mathrm{\Delta }𝐓^2`$. For $`\mathrm{\Delta }=0`$, we have $`\mathrm{\Delta }𝐓^2=T`$. For $`\mathrm{\Delta }>0`$, the diffuseness of the nucleon numbers makes a contribution to $`\mathrm{\Delta }𝐓^2`$ which is largest for $`T=0`$. $`\mathrm{\Delta }𝐓^2`$ is an even function of $`T`$ and for $`\mathrm{\Delta }>0`$ analytic at $`T=0`$. This results in the hyperbola-like behaviour seen in the inset of Fig. 3 of Ref.  and reflected in Fig. 2. For $`\mathrm{\Delta }=0`$, the first term in $`E_2^{}`$ has a $`T`$-dependent part $`GT/2`$, which arises from blocking of the neutron-proton pairing force when nucleons are promoted from proton to neutron levels. For $`\mathrm{\Delta }>0`$, the asymptotic slope of the hyperbola gets a corresponding contribution from this term. Since $`G=.4\text{MeV}\eta +\kappa =3.3\text{MeV}`$, this is small compared to the contribution of the second term. Evidently, this discussion of $`E_2^{}`$ applies analogously to the true Fock term without the term proportional to $`\eta `$.

As seen from Fig. 1, the RPA contribution to the symmetry energy is in the pairing plus symmetry force model largely taken into account by adding $`\mu /2`$ to the self-consistent energy $`E_0`$. This corresponds to the approximation employed by Frauendorf and Sheikh, who added $`\mu /2`$ to an isocranked Hartree-Fock-Bogolyubov energy . As pointed out by these authors, doing so is essentially equivalent to assigning the energy calculated for $`T_\alpha =T+1/2`$ or $`\sqrt{T(T+1)}`$ to a state with isospin $`T`$, where $`T_\alpha `$ is some isospin coordinate. In fact, if $`E_0(T_\alpha )`$ is the calculated energy, we have $`E_0\left(\sqrt{T(T+1)}\right)E_0(T+1/2)E_0(T)+E_0^{}(T)/2=E_0(T)+\mu /2`$. Such procedures are common in ordinary cranking calculations. See references in Ref. . In a context of isocranking, Satula and Wyss applied the constraint $`T_x=\sqrt{T(T+1)}`$ . In the pairing plus symmetry force model, the $`T`$-dependence of the RPA correlation energy is seen to be taken into account by such recipes only in a crude approximation: $`E_2`$ has a $`T`$-dependence additional to that of $`\mu /2`$. It is in fact this additional $`T`$-dependence which makes $`x>1`$.

The presence of a strong isovector pair field is essential for getting $`x1`$. For $`G=0`$ the Hartree-Fock approximation is in fact exact, and we have $`EE_{T=0}=(\eta T^2+\kappa T(T+1))/2=(\eta +\kappa )/2\times T(T+x)`$, where $`x=\kappa /(\eta +\kappa )=.7`$ for $`\eta =1.1\text{MeV}`$. The leading contribution of the pairing force in perturbation theory $`G𝐏^{}𝐏=G(T3\mathrm{\Omega })/2`$ adds a linear term $`GT/2`$ to the symmetry energy, but as long as $`G`$ is below the threshold of pair condensation $`G\eta /(2+\mathrm{log}(\mathrm{\Omega }/2))`$ it is much smaller than the term $`\eta T/2`$ required to get $`x=1`$. A pair gap $`\mathrm{\Delta }`$ which is comparable to the single-nucleon level spacing $`\eta `$ is thus necessary in order to make the quasiparticle energies and the RPA frequencies except $`\mu `$ so independent of $`T`$ that $`\mu /2`$ becomes the dominant term in $`E_2E_{2,T=0}`$.

In conclusion, the pairing plus symmetry force model with equidistant single-nucleon levels introduced in Ref.  gives in the RPA for $`A=48`$ a symmetry energy proportional to $`T(T+x)`$, where $`x1.025`$. The large isovector pair field, whose strength is inferred from the observed odd-even mass staggering, is essential for this result. The RPA contribution to the symmetry energy is roughly equal to $`\mu /2`$, where $`\mu `$ is the isocranking angular velocity, but there are significant corrections to this approximation. It is in fact these corrections which make $`x>1`$. The RPA contribution to the symmetry energy cannot be simulated by adding to the Hartree-Fock-Bogolyubov energy a term proportional to the isospin variance in the Bogolyubov quasiparticle vacuum.


text:
# An empirical test of the theoretical population corrections to the Red Clump absolute magnitude

## 1 Introduction

In recent years much work has been devoted to studying the suitability of Red Clump (RC) stars as a distance indicator. These core helium-burning stars are essentially in an identical evolutionary phase to those which make up the horizontal branch in globular clusters. However, in intermediate and higher metallicity systems only the red end of the distribution is seen, forming a clump of stars in the colour-magniutude diagram - hence the name. In an infrared study of metal-rich globular clusters, Kuchinski et al. (1995) realised that the absolute $`K`$-band magnitude of the horizontal branch (effectively, the RC) displayed very little cluster-to-cluster variation, and therefore suggested it could potentially be used as a standard candle for distance determinations. However, at that time, no absolute calibration of the RC magnitude was available and so distances could only be measured differentially, with respect to a specific cluster. Interest in the RC as a potential extra-galactic distance indicator was rekindled a few years later due to the availability of accurate $`Hipparcos`$ (ESA 1997) parallaxes (and hence absolute magnitudes) for several hundred local RC stars, enabling a very precise calibration of their mean absolute magnitude. In fact, RC stars are currently the only extra-galactic distance indicator which can be tied precisely to the $`Hipparcos`$ distance scale, due to the lack of a significant number of RR Lyraes with accurate parallax determinations (see e.g. Fernley et al. 1998) and the current uncertainties in the metallicity dependence of the Cepheid period–luminosity relation (Feast 2003, and references therein).

Paczyński & Stanek (1998) made the first attempt to provide a precise measurement of the local RC absolute magnitude from $`Hipparcos`$ data by fitting the $`I`$-band magnitudes of RC stars to a Gaussian function, superimposed on a background distribution of red giant branch (RGB) stars, to find a peak value of $`M_I^{RC}=0.28`$. This early work focussed on the $`I`$-band since, observationally, the $`I`$-band magnitudes of RC stars appeared to show no dependence on colour, unlike those in the $`V`$-band, hence it was thought that the $`I`$-band magnitude of RC stars may provide a truly ‘standard’ candle. The initial value from Paczyński & Stanek (1998) was subsequently revised slightly to $`M_I^{RC}=0.23`$ by Stanek & Garnavich (1998). This was then used by Stanek, Zaritsky & Harris (1998) to make the first direct comparison with the RC in the LMC, yielding a dereddened distance modulus of $`\mu _{0,LMC}=18.065\pm 0.031_{\mathrm{statistical}}\pm 0.09_{\mathrm{systematic}}`$. Contemporaneously, an almost identical result of $`\mu _{0,LMC}=18.08\pm 0.03_{\mathrm{statistical}}\pm 0.12_{\mathrm{systematic}}`$ was found by Udalski et al. (1998) using the same method. The implied LMC distance from these results is short, even by the standards of the ‘short’ distance scale, and they are in significant disagreement with the LMC distance modulus of $`\mu _{0,LMC}=18.50\pm 0.10`$ adopted by the $`HubbleSpaceTelescope`$ Extragalactic Distance Scale Key Project (Freedman et al. 1994).

It was soon pointed out that models of core helium-burning stars predict that the RC luminosity is dependent on both age and metallicity, causing differences of up to 0.6 mag in the mean value of $`M_I^{RC}`$ between different stellar populations (Cole 1998; Girardi et al. 1998). The implication therefore, is that any distance derivation to an extragalactic system should first account for population differences between the local RC and the RC of the system in question, for the results to be meaningful. In an attempt to address this problem, Udalski (2000) examined the metallicity dependence of $`M_I^{RC}`$ in the solar neighbourhood from the empirical data. He found a weak dependence of $`M_I^{RC}`$ on \[Fe/H\] at the level of $``$ 0.13 mag $`\mathrm{dex}^1`$, and consequently derived an $`I`$-band correction of $``$0.07 mag for the LMC. Subsequent work has focussed on tightening the constraints on the metallicity dependence of $`M^{RC}`$ in several pass-bands, whilst at the same time, investigating the effects of age.

Based on the models of Girardi et al. (2000), predicted mean RC magnitudes for a large range of ages (0.5$``$12 Gyr) and metallicities ($`1.7[\mathrm{Fe}/\mathrm{H}]+0.2`$) are given by Girardi & Salaris (2001, hereafter GS01, $`V`$ and $`I`$-band) and Salaris & Girardi (2002, hereafter SG02, $`K`$-band). In order to test these predictions the effects of age and metallicity must be examined separately. To this end, a number of recent studies have examined small samples of Galactic open clusters, each of single age and chemical composition, enabling a comparison of the observed RC magnitudes with those predicted from theory (Sarajedini 1999; Twarog, Anthony-Twarog & Bricker 1999; GS01; Grocholski & Sarajedini 2002). Qualitatively, these studies find broad agreement with the models in all 3 pass-bands ($`V`$, $`I`$ and $`K`$), confirming that ‘population corrections’ to the RC absolute magnitude are necessary before it can be used as an accurate extragalactic distance indicator. It is important to note here that the metallicity dependence found by Udalski (2000) is derived from local stars in the disc and cannot be regarded as universally applicable since, if we accept the need for population corrections at whatever level, the star formation history and age-metallicity relation are likely to be key factors in determining the metallicity dependence of the RC magnitude. In fact, GS01 model the RC in the solar neighbourhood and reproduce very well the metallicity dependence found by Udalski, giving further weight to the applicability of the models.

In this work, we give a more quantitative assessment of the accuracy of the theoretical predictions for the RC population corrections, over a range of ages and metallicities. We also aim to improve on previous work in several ways. Firstly, determination of the RC absolute magnitude for a particular cluster requires knowledge of the cluster distance. All previous studies involving samples of galactic open clusters have used distances which have been derived using methods based on fitting the cluster sequences to theoretical isochrones. For example, Sarajedini (1999) employs cluster distances which are derived relative to that of M67, the distance to M67 having been determined by fitting to theoretical isochrones. These methods can potentially introduce systematic errors which depend on the particular choice of isochrones. In this work, we employ the purely empirical MS-fitting method, based on a large sample of local field dwarfs with accurate $`Hipparcos`$ parallaxes, developed in Percival, Salaris & Kilkenny (2003, hereafter Paper I), thus removing a potentially large source of systematic error. Furthermore, we ensure complete consistency in the analysis by studying $`V`$-, $`I`$\- and $`K`$-band data simultaneously, utilizing the same distance modulus, age estimate and reddening for each cluster in all 3 bands.

Secondly, using Monte Carlo methods, we quantify the effects of sampling on the luminosity of the RC in clusters of different ages and metallicities, in terms of its (empirically) measured magnitude and the likely error. The same simulations also allow us to quantify any offsets between mode, median and mean values for the RC absolute magnitude, $`M_\lambda ^{RC}`$. It should be noted that the model values tabulated in GS01 and SG02 are the predicted mean values for $`M_\lambda ^{RC}`$ whilst, from empirical data, generally the median of the distribution is measured. These values are expected to differ since the distribution of magnitudes in the RC is not symmetric – the mode of the magnitude distribution is centred on the zero age horizontal branch, and the distribution has a tail towards brighter magnitudes as the RC stars evolve to higher luminosities.

A third important element in this work, particularly when considering the population effects in the LMC, is the inclusion of 4 clusters with ages $``$ 2 Gyr. These young clusters are important since, for populations with recent star formation (like the LMC), a large fraction of RC stars have ages in the range 1$``$2 Gyr. The models predict large, and rapid, changes in $`M^{RC}`$ in all 3 pass-bands for stars with ages $`<`$ 2 Gyr and so it is vitally important to test these predictions empirically, and quantitatively assess their validity.

The layout of the paper is as follows: In Section 2 we review the basic method of using the RC as a distance indicator, in Sect. 3 we give details of the clusters used in this study and in Sect. 4 we describe the methods involved in deriving cluster distances, and determining the RC magnitude, both empirically and from the theoretical models. Sect. 5 presents our results and in Sect. 6 we discuss some implications of these results and make our conclusions.

## 2 The RC as a Distance Indicator – Basic Method

In order to utilize the Red Clump as an extragalactic distance indicator, once the apparent magnitude of the RC, $`m_\lambda ^{RC}`$, is measured, the distance modulus, $`\mu _0=(mM)_0`$, to a particular galaxy is given by:

$$\mu _{0,galaxy}=m_\lambda ^{RC}M_\lambda ^{RC}A_\lambda +\mathrm{\Delta }M_\lambda ^{RC}$$
(1)
where $`m_\lambda ^{RC}`$ is the apparent magnitude of the RC in a particular pass-band, $`M_\lambda ^{RC}`$ is the ‘zero-point’ absolute RC magnitude given by the local RC stars, $`A_\lambda `$ is the extinction, and $`\mathrm{\Delta }M_\lambda ^{RC}`$ is the population correction to the RC magnitude.

The most recent determinations of the local RC absolute magnitude, $`M_\lambda ^{RC}`$ are given by Alves et al. (2002), and are $`0.73\pm 0.03`$, $`0.26\pm 0.03`$ and $`1.6\pm 0.03`$ mag in, respectively, the $`V`$, $`I`$ and $`K`$ pass-bands. $`\mathrm{\Delta }M_\lambda ^{RC}`$, the population correction, is defined as the difference between the local RC absolute magnitude and that of the RC of the system in question i.e. $`\mathrm{\Delta }M_\lambda ^{RC}M_\lambda ^{RC}(local)M_\lambda ^{RC}(galaxy)`$. As we have already seen, theory predicts that $`\mathrm{\Delta }M_\lambda ^{RC}`$ is dependent on both age and metallicity.

To test for and quantify these dependencies, we need to employ a sample of galactic open clusters i.e. simple populations, each of single age and single metallicity, so that the effects of age and metallicity can be disentangled. From eq. 1, it can be seen that the population correction, $`\mathrm{\Delta }M_\lambda ^{RC}`$, for a particular cluster will be given by:

$$\mathrm{\Delta }M_\lambda ^{RC}=\mu _{0,cluster}m_\lambda ^{RC}+M_\lambda ^{RC}+A_\lambda $$
(2)
In this equation, $`m_\lambda ^{RC}`$, $`M_\lambda ^{RC}`$ and $`A_\lambda `$ are known quantities from the literature (see sect. 3 for full details) and the distance moduli to the individual clusters, $`\mu _{0,cluster}`$ are to be determined.

As already noted, all previous studies have relied on cluster distances derived from MS fitting to theoretical isochrones (i.e. using stellar models), which potentially introduces systematic errors. The aim of this work is to quantify the population corrections via a completely model independent method by utilizing our purely empirical MS-fitting procedure (Paper I) to derive the cluster distances. This method employs a sample of 54 local unevolved field dwarfs with accurate $`Hipparcos`$ parallaxes (errors are generally $`<5\%`$) for which we acquired new $`BVI_C`$ photometry. From fitting to this sample in both the $`(BV)`$ and $`(VI)`$ colour planes we derive a Hyades distance modulus of $`(mM)_0=3.33\pm 0.06`$ – exactly reproducing the distance directly measured from parallaxes of the Hyades stars (Perryman et al. 1998). Hence we are able to tie the absolute magnitude of the RC in each of the open clusters in our sample to the $`Hipparcos`$ distance scale. Since the absolute magnitude of the local RC is also determined from $`Hipparcos`$ parallaxes, systematic errors are minimized and we can provide a completely self-consistent assessment of the age and metallicity dependencies for the RC population corrections.

## 3 Choosing the Cluster Sample

The two most important criteria in choosing open clusters for inclusion in this work are that, firstly, each cluster must have a very well defined main sequence in the HR-diagram so that ridge lines can be derived and, secondly, the photometry (which is all existing in the literature) must go deep enough to apply our empirical MS-fitting method. This method uses only the unevolved lower portion of the MS in the fitting procedure, to ensure that there are no evolutionary effects – hence the cluster photometry must extend down to at least $`M_V=7.0`$. Since we have $`B`$, $`V`$ and $`I_C`$ photometry for our field dwarf sample used in the MS fitting, we can utilize clusters with either $`(BV)`$ or $`(VI)`$ photometry or, indeed, both.

We took, as a starting point, the sample of 8 open clusters from Sarajedini (1999, hereafter S99), however we discounted NGC 6791 and NGC 6819 as their main sequence photometry was not deep enough to apply our MS-fitting method (see Section 4.1 for more details). Thus we included M67, NGC 188, NGC 7789, Berkeley 39, NGC 2204 and Mel 66 – the sources of photometry for these 6 clusters are as listed in S99, their Table 1. We added 2 more young clusters for which there is deep photometry available in the literature, namely NGC 2204 (Kassis et al. 1997, $`(BV)`$ and $`(VI)`$ data) and NGC 2420 (Anthony-Twarog et al. 1990, $`(BV)`$ data only).

Cluster abundances (and their quoted errors) were taken from the catalogue of Gratton (2000), ensuring that the metallicities for the sample clusters and for the field dwarfs used in the MS-fitting are on an homogeneous scale (see discussion in Paper I). We took cluster reddenings from Twarog, Ashman & Anthony-Twarog (1997), assuming a 0.02 mag error in $`E(BV)`$ for all. Absolute ages (and associated errors) were estimated from comparison of the dereddened cluster main lines, and their turn off absolute magnitudes, with Girardi et al. (2000) isochrones - the absolute magnitude of the turn off having been determined using the distance moduli derived in this paper. We note here that our derived ages do not differ significantly from those found by S99, who used the Bertelli et al. (1994) isochrones to estimate cluster ages.

These 8 open clusters span a metallicity range of $`0.44[\mathrm{Fe}/\mathrm{H}]+0.02`$ and have ages from 1$``$7.5 Gyr. We also included in the sample the metal-rich Galactic Globular Cluster 47 Tuc ($`[\mathrm{Fe}/\mathrm{H}]=0.7`$), since we have an empirically determined MS-fitting distance (and hence an age estimate) for this cluster which is also tied to the $`Hipparcos`$ distance scale (Percival et al. 2002). Photometry for 47 Tuc is from Stetson’s ‘secondary’ standards available through the Canadian Astronomy Data Centre web site (http://cadcwww.hia.nrc.ca/standards/ – see Stetson 2000 for details). This extends the age range of the cluster sample up to 11 Gyr. In addition, 6 of these clusters have $`K`$-band data from the Two Micron All Sky Survey (2MASS) Point Source Catalogue (http://irsa.ipac.caltech.edu/) which is used by Grocholski & Sarajedini (2002) to determine the absolute $`K`$-band magnitude of the RC – these are; NGC 2204, Berkeley 39, NGC 2420, NGC 2477, M67 and 47 Tuc. The basic cluster parameters are listed in Table 1.

## 4 Method

### 4.1 Empirical MS fitting

Before the MS-fitting method can be applied to derive the cluster distance, the MS ridge line must be determined. The procedure we adopted consisted of making 0.2 mag cuts in $`V`$ across the MS and then plotting histograms in colour ($`(BV)`$ or $`(VI)`$) to find the peak value. This method clearly identifies the ridge line for a cluster and also clearly shows the location of a binary sequence, if present. The resultant points are fitted to a polynomial (usually cubic) function.

Cluster distances are derived using our empirical MS-fitting method, which we briefly outline here. The method is fully described in Paper I, to which we refer the interested reader for more details. In Paper I we presented new $`BVI`$ (Johnson-Cousins) photoelectric photometry for a sample of 54 field dwarfs with accurate $`Hipparcos`$ parallaxes, which span a metallicity range of $`0.4[\mathrm{Fe}/\mathrm{H}]+0.3`$. Absolute magnitudes for these stars are in the range $`5.4M_V7.0`$, thus they are unevolved zero-age MS stars and their location in a Colour-Magnitude diagram is totally insensitive to any age differences between the clusters which we are fitting to and field dwarf sample.

The basic MS-fitting method consists of constructing a template MS from the field dwarf sample by applying colour shifts to the individual stars to account for the differences in metallicity between the field stars and the cluster. This template is then shifted in $`V`$ to match the dereddened and extinction-corrected cluster main line, the difference in magnitude being equal to the absolute distance modulus $`(mM)_0`$.

The procedure used to calculate the magnitude of the metallicity dependent colour shifts for the field dwarfs relies on first establishing that the shape of the MS is insensitive to \[Fe/H\] in the narrow range of metallicities and magnitudes we are dealing with. Next, we determine the colour that each field dwarf would have at a fixed magnitude of $`M_V=6`$, which we call $`(BV)_{M_{V=6}}`$ (or $`(VI)_{M_{V=6}}`$), using the slope of the Hyades MS as a reference slope. Finally, we calculate the derivative $`\delta (BV)_{M_{V=6}}/\delta [Fe/H]`$ (or $`\delta (VI)_{M_{V=6}}/\delta [Fe/H]`$). Because of the constant shape of the MS across the metallicity range, this derivative is appropriate for the whole magnitude range spanned by the field dwarfs, and hence colour shifts are applied to each star in the sample, at their observed magnitudes, to construct a template MS at the metallicity of the cluster in question.

When fitting to the cluster main lines to derive the cluster distance, weighted errors are used for the field dwarfs which include both photometry errors and magnitude errors due to error on the parallax (where $`\mathrm{\Delta }M_V=2.17(\mathrm{\Delta }\pi /\pi )`$). We also account for errors on the metallicity of the both the field dwarfs and the cluster, and the best-fitting distance is found by minimizing $`\chi ^2`$. For clusters with photometry in both colour indices, dereddened distances are taken to be an average from the $`(BV)`$ and $`(VI)`$ fits, which we note here are always consistent with each other, within their respective errors.

The fitting procedure is illustrated in Figure 1 which shows the $`V/(BV)`$ CMD for NGC 2477 (data from Kassis et al. 1997). The left panel shows the dereddened and extinction corrected cluster photometry overlaid with the derived ridge line for the MS. The RC stars are contained within the box spanning $`11.1V_012.1`$ and $`0.85(BV)_01.1`$. The panel on the right shows the cluster ridge line, shifted in $`V`$-magnitude to match the template MS created from the field dwarf sample, hence yielding the best-fitting distance modulus. The data points are the template MS, which has been created from the field dwarf sample by ‘correcting’ their colours to match the metallicity of the cluster ($`[\mathrm{Fe}/\mathrm{H}]=0.0`$ in this case), as described above. The dereddened distance moduli determined by this method for all the clusters in the sample are listed in Table 1.

### 4.2 RC magnitudes - observed

Apparent RC magnitudes, $`m_\lambda ^{RC}`$, are taken as listed in S99 ($`V`$\- and $`I`$-bands only) for M67, NGC 188, NGC 7789, Berkeley 39, NGC 2204 and Mel 66. For NGC 2204 and NGC 2420 we calculated the RC magnitude using exactly the same method employed by S99, i.e. finding the median magnitude of the stars contained in a box centered on the RC. As noted by S99, this median value is virtually insensitive to the exact location of the box, a fact confirmed by us from Monte Carlo simulations of the RC (see Sect. 4.3). $`K`$-band RC magnitudes, determined using the same method, are taken from Grocholski & Sarajedini 2002 (note that these are listed as extinction corrected absolute magnitudes, $`M_K^{RC}`$, so we have used their quoted reddening and distance modulus values to calculate the apparent $`K`$ magnitudes, $`m_K^{RC}`$).

Thus, for each cluster the absolute RC magnitude, $`M_\lambda ^{RC}`$ is given by:

$$M_\lambda ^{RC}(cluster)=m_\lambda ^{RC}A_\lambda (mM)_0$$
(3)
where $`m_\lambda ^{RC}`$ is the apparent RC magnitude as described above, $`A_\lambda `$ is the extinction and $`(mM)_0`$ is the dereddened distance modulus as found in Sect. 4.1.

The population correction, $`\mathrm{\Delta }M_\lambda ^{RC}`$, is then given by

$$\mathrm{\Delta }M_\lambda ^{RC}=M_\lambda ^{RC}(local)M_\lambda ^{RC}(cluster)$$
(4)
where $`M_\lambda ^{RC}(local)`$ are the values given by Alves (2002), for $`V`$, $`I`$ and $`K`$ (see Sect. 2).

### 4.3 RC magnitudes - theoretical

Mean RC magnitudes generated from the Girardi et al. (2000) models are tabulated in GS01 ($`V`$ and $`I`$) and SG02 ($`K`$) for a range of ages and metallicities (see Sect. 1). These effectively correspond to cluster RC magnitudes, $`M_\lambda ^{RC}(cluster)`$, since each tabulated value represents the RC magnitude for a single age and metallicity.

For the mean magnitude of the local RC, $`M_\lambda ^{RC}(local)`$, we use the values calculated by GS01 and SG02. GS01 and SG02 employ a complete population synthesis algorithm which, using the stellar models of Girardi et al. (2000), produces a synthetic CMD, from which the RC luminosity function can be extracted. The solar neighbourhood was modelled using the Star Formation Rate (SFR) and Age Metallicity Relation (AMR) from Rocha-Pinto et al. (2000a,b). We will discuss the effects of the particular choice of SFR and AMR in Section 6. The model population correction is then calculated as for the empirical data i.e. $`\mathrm{\Delta }M_\lambda ^{RC}=M_\lambda ^{RC}(local)M_\lambda ^{RC}(cluster)`$. It is important to note at this point that we are not relying on the *absolute* values for the RC magnitudes from the models since the population correction is a differential quantity, i.e. it is the *difference* between the local RC magnitude and the RC magnitude of the population under scrutiny.

The tabulated model values give the mean properties of the RC stars – i.e. the mean magnitude of all stars in the simulation which are still undergoing core helium burning. RC stars start their lives on the ZAHB, producing a peak in their luminosity function, and then evolve on shorter timescales to brighter magnitudes. Consequently, the mode of the magnitude distribution is offset from the mean value, in the sense that the mean magnitude is brighter. Obviously, from the simulations, all the RC stars are accounted for and any value can be measured – mean, median or mode. For real data sets the RC is often sparsely populated, due to poor sampling, and the empirical median values derived as described in Sect. 4.2 reflect the mode rather than the mean magnitude.

We tested for the offset between mean and median magnitudes by first simulating RCs at ages from 1–10 Gyr in $`V`$, $`I`$ and $`K`$. This was done by populating the Girardi et al. (2000) isochrones for solar metallicity using a Salpeter IMF, and then taking all stars still in their core helium-burning phase. These stars are all within 1.0 mag of the ZAHB and, in the particular simulations used, we generate a total parent sample of 4000–5000 RC stars. (We note here that the particular choice of IMF does not influence the results at all, since the evolutionary time-scale for these stars is very short, and so the mass is nearly constant in the region of interest). We then randomly extracted subsamples of 5 and 50 stars from the parent sample and calculated the median magnitude in all 3 pass-bands for each subsample. This process was repeated 70 times and the resultant median magnitudes plotted as histograms. The peak of each distribution (i.e. for the subsamples of 5 and the subsamples of 50) yields the most likely measured median value for a real data set, and the dispersion (i.e. 1-$`\sigma `$) of the distribution represents the error associated with the (small) sample size. The median values were found to be identical (within 0.01 mag) for the simulated subsamples of 5 and 50, and whilst the errors are obviously larger for the smaller sample size, they are surprisingly stable – errors are generally no more than 0.05 mag for a sample of 5, except for the 1 Gyr simulation which produces a 1-$`\sigma `$ error of $``$ 0.1 mag. However, this is not a cause for concern here, as all the young clusters used in this study have a well populated RC, which produces much smaller errors.

The mean–median offsets are 0.04–0.05 mag for ages $`>`$ 2 Gyr in all 3 filters, whilst at 1 Gyr, the offset is 0.07–0.08 mag. Our empirical RC absolute magnitudes, $`M_\lambda ^{RC}`$, were therefore adjusted using the appropriate mean$``$median offset before calculating the population corrections, $`\mathrm{\Delta }_\lambda ^{RC}`$, and comparing with the theoretical values. The empirical values of $`M_\lambda ^{RC}`$ are listed in Table 1 – the quoted errors include the error on the cluster distance modulus, the statistical error on the measured RC magnitude, described above, and the error induced by a 0.02 mag uncertainty on the cluster reddening, added in quadrature.

## 5 Results

Fig. 2 shows the empirically derived population corrections, $`\mathrm{\Delta }M_\lambda ^{RC}`$, in $`V`$ (9 pts), $`I`$ (8 pts) and $`K`$ (6 pts) overplotted with the theoretical models of Girardi et al. (2000) for \[Fe/H\] between $``$0.7 and $`+`$0.2. The data points are our empirically determined values and the vertical error bars include the error on our distance moduli (which includes error on cluster metallicity), the statistical error on the RC magnitude (detailed in Sect. 4.3), and the error on the RC magnitude induced by uncertainty in the cluster reddening (0.02 mag), all added in quadrature. The most important feature of this plot is the fact that the empirically determined population corrections follow the trend of the theoretical models in all 3 pass-bands, particularly at ages $`<`$ 2 Gyr, where the models predict a strong dependence on age.

We stress again that we are not relying on *absolute* values from the models, which may differ from the empirical ones (e.g. for local RC magnitude, see GS01). The models are always being used differentially, remembering that the population correction, $`\mathrm{\Delta }M_\lambda ^{RC}`$, is the *difference* between the local RC magnitude and the RC magnitude of the cluster in question.

In order to quantify the goodness-of-fit between the models and the empirical data, we calculated the residuals of $`\mathrm{\Delta }M_\lambda ^{RC}`$ in the sense $`modelobserved`$. These residuals are plotted in Figs. 3 and 4 against, respectively, age and metallicity. The vertical error bars (i.e. in magnitude only) are as in Fig. 2, but have been amended to account for the error on the cluster age. It can be seen from Fig. 2 that the error on the cluster age will cause an asymmetric error on $`\mathrm{\Delta }M_\lambda ^{RC}`$, particularly for ages $`<`$ 2 Gyr, and this is reflected in the error bars in Figs. 3 and 4. These plots clearly show that the residuals in $`\mathrm{\Delta }M_\lambda ^{RC}`$ display no statistically significant trends with either age or metallicity in any of the 3 pass-bands.

The most discrepant point in the $`V`$\- and $`I`$-bands (i.e. with the largest residuals) is Melotte 66. This is the only cluster for which there is no $`(BV)`$ data for the main sequence and so the distance modulus has been derived from the $`(VI)`$ data alone. We note, however, that there is a small amount of photoelectric $`UBVRI`$ photometry for this cluster from Twarog, Anthony-Twarog & Hawarden (1995). Comparison of the photoelectric data with the CCD data used in this study (from Kassis et al. 1997) indicates that there may be small zero-point offsets in the CCD photometry which would have the effect of making the $`(VI)`$ colours too red by 0.04–0.05 mag. If this were the case, the derived distance would be too short by $``$ 0.2 mag, causing a similar discrepancy in the RC magnitude, and therefore accounting for the discrepancy between the model and observed values for the population correction. Removing this cluster from our analysis, the mean residuals in $`\mathrm{\Delta }M_\lambda ^{RC}`$ are; +0.03 $`\pm `$0.07 in $`V`$, +0.07 $`\pm `$0.06 in $`I`$, and 0.00 $`\pm `$0.07 in $`K`$, where the quoted errors are the 1-$`\sigma `$ errors. Hence we find no statistically significant offsets in $`V`$ and $`K`$ and a very marginal offset in $`I`$, significant only at the 1.2-$`\sigma `$ level.

We note that if we included Melotte 66 in this analysis the mean residuals in $`V`$ and $`I`$ would become +0.05 $`\pm `$0.09 and +0.09 $`\pm `$0.06 respectively (those in $`K`$ would remain unchanged) and hence our conclusions would essentially remain unaltered.

## 6 Discussion and conclusions

A vital factor in determining the RC population corrections is the absolute magnitude of the local RC, for which we have very well determined observational values. Deriving the same quantities from the models requires knowledge of the local SFR and AMR, which are used as input parameters in the population synthesis algorithm (see GS01 for details). We checked the effects of assuming a different SFR and/or AMR by modelling the local RC using the AMRs of Carlberg et al. (1985) and of Edvardsson et al. (1993). We also tested the effect of using a constant star formation rate. In all cases, the model values for the mean RC magnitude, $`M_\lambda ^{RC}`$, were indistinguishable from those found when assuming the Rocha-Pinto et al. (2000a,b) relationships. This is largely because the local RC is dominated by stars with ages $`<`$ 2 Gyr (see in depth discussion in GS01) and the different AMRs all predict similar metallicities at these young ages. Hence the precise details of the star formation history do not significantly affect the local RC magnitude.

The local RC simulated in GS01 (using the Rocha-Pinto 2000 SFR and AMR) yields a straight mean metallicity of $`[\mathrm{Fe}/\mathrm{H}]=0.04`$, whilst the observed RC stars in the solar neighbourhood appear to have a mean \[Fe/H\] of $``$ -0.18 (see, e.g. Alves 2000). Metallicities used by Alves for the local RC stars are from McWilliam (1990) – the same ones used by Udalski (2000) to derive a metallicity dependence for $`M_I^{RC}`$. This offset between the model and observed mean metallicity for the RC stars has already been noticed and discussed in GS01. McWilliam himself notes that his metallicities for the RC stars may be systematically too low, especially for the highest metallicities, due to the (then) lack of suitable model atmospheres for super-solar metallicities. Luck & Challener (1995) determined abundances for 55 G and K field giants from high resolution spectra – 45 of these stars being in common with the McWilliam (1990) sample. Luck & Challener find a constant offset of 0.12 dex between their abundances and those of McWilliam, in the sense that the McWilliam (1990) values are lower. From the most recent catalogue of \[Fe/H\] determinations for FGK stars by Cayrel de Strobel, Soubiran & Ralite (2001), we made a comparison of all the stars in the McWilliam (1990) sample which had abundance determinations by other authors – this yielded a sample of 89 stars with 152 measurements from different sources. The average difference between the McWilliam abundances and all the others in the literature is 0.13 dex, with a dispersion of 0.15 dex. However it is evident that the McWilliam metallicities are increasingly discrepant at the highest metallicities. Splitting the sample into 2 groups with $`[\mathrm{Fe}/\mathrm{H}]_{\mathrm{other}}<0.2`$ and $`[\mathrm{Fe}/\mathrm{H}]_{\mathrm{other}}>0.2`$ and comparing with the McWilliam values reveals an offset of 0.10 dex for the lower metallicity group and an offset of 0.31 dex for the higher metallicity group, both with a dispersion of 0.14 dex, in the sense that the McWilliam values are lower. Hence, the apparent discrepancy between the average metallicity of the local RC predicted from the models and that derived from observational data, appears to be a problem with the zero-point of the metallicity scale derived by McWilliam (1990).

The results of our analysis fully confirm the theoretical population corrections obtained by GS01 and SG02 for single-age, single-metallicity stellar populations, which are the building blocks for computing corrections to composite systems. Remarkably, the minimum of the RC luminosity for ages around 1 Gyr predicted by theory, corresponding to the transition from electron degenerate to non-degenerate helium cores, is also clearly observed. This is a further positive test for the accuracy of the physics employed in the stellar models. A larger sample of clusters with accurate multicolour photometry for the lower MS could confirm whether the small offset in the $`I`$-band is real. For use as an accurate distance indicator, $`K`$ is the best single pass-band to use as the effects of reddening are minimised and, provided ages are less than $``$ 7 Gyr, the metallicity dependence is very small.

It is evident from Fig. 2 that, in general, there is no simple relationship between the population correction in a given passband and metallicity or age. Also, it is only in the $`K`$-band, and specifically for metallicities close to solar and ages lower than a few Gyr, that the population corrections are very small (see SG02 for a discussion on this point). This means that there is no way, a priori, to determine generic empirical population corrections that are a function of only the metallicity of the population under scrutiny (this point has already been discussed in depth by GS01). In the case of composite stellar populations it is necessary to have a good determination of their Star Formation Rate (SFR) and Age Metallicity Relationship (AMR) and apply population synthesis algorithms in order to obtain the appropriate corrections (alternatively, the corrections can be derived using the simpler analytical method described in GS01). The accuracy of the derived corrections will depend entirely on the accuracy of the determined SFR and AMR, since they are combinations of the single-age, single-metallicity corrections, which we have now proved to be accurate. For systems where no precise AMR and SFR are available, the use of empirical calibrations – which are necessarily of limited validity since they are dependent on the particular properties of the calibrating sample – can be very dangerous because it can lead to errors of several tenths of a magnitude on the derived distance moduli.

The very limited range of validity of empirical corrections based only on the stellar metallicity is illustrated by the valuable data published by Pietrzyński, Gieren & Udalski (2003, hereafter PGU03). These authors present new near-infrared $`J`$\- and $`K_\mathrm{S}`$-band data for several fields in the LMC, SMC, and the Fornax and Carina dwarf galaxies and, in addition, they utilize $`I`$-band data in the literature for the same systems. PGU03 are therefore able to determine accurately the $`K`$\- and $`I`$-band brightness of the RC in these 4 systems and they go on to derive empirical corrections to the RC brightness, based solely on differences in mean metallicity. After applying these corrections (which are normalized to the LMC mean metallicity) to the RCs in the SMC, Carina, Fornax and the solar neighbourhood, PGU03 look at the difference between the $`I`$ and $`K`$ brightness of the RC for each system. If the corrections are adequate, one expects that this quantity will be constant for all the systems, however PGU03 find a discrepancy of 0.23 mag for the solar neighbourhood value, which they ascribe to zero-point errors in the photometry of the local RC stars.

This difference between the $`I`$\- and $`K`$-band absolute magnitude of the RC is an ideal parameter on which we can test the theoretical models further. The predicted metallicity dependencies for the population correction, $`\mathrm{\Delta }M_\lambda ^{RC}`$, are very different between the 2 passbands, that in $`I`$ being much stronger than in $`K`$. Furthermore, the dependencies work in opposite directions – for increasing metallicity, $`\mathrm{\Delta }M_I^{RC}`$ becomes smaller, whilst $`\mathrm{\Delta }M_K^{RC}`$ becomes larger – and are not constant with increasing age. Therefore, as an experiment, we examined the difference between the $`I`$-band and $`K`$-band magnitude of the RC predicted from theory for the LMC, SMC, Carina and the solar neighbourhood. It is important to realise that in case of the theoretical population corrections, consistency with the local RC is automatically achieved, since the corrections are computed relative to the solar neighbourhood values. We have adopted the values for the LMC and SMC provided by SG02; for Carina, we have simulated the RC using the same population synthesis algorithm, assuming the SFR used in GS01 and adopting the most recent determination of the AMR from Tolstoy et al. (2003). We did not consider the case of Fornax simply because existing determinations of the SFR are highly uncertain, thus preventing an accurate assessment of the appropriate theoretical population corrections.

Figure 5 shows a comparison of the observed difference between the $`I`$\- and $`K`$-band brightness of the RC, for the LMC, SMC, Carina and the solar neighbourhood, and the corresponding values predicted from theory. (We note here that the $`K_\mathrm{S}`$-band magnitudes from PGU03 have been transformed to the Bessell & Brett system, which are those used in the theoretical models). The result of the comparison is excellent, with a good agreement, within the errors, between theory and observations. This confirms the accuracy of the theoretical corrections, which simultaneously reproduce the LMC, SMC and Carina data and satisfy the solar neighbourhood values without the need to invoke photometric errors of the order of 0.20 mag. As a consequence, with the theoretical population corrections one obtains full consistency between the RC $`I`$\- and $`K`$-band distances to the LMC derived from the empirical values tabulated by PGU03; in fact, by using the local RC calibration by Alves et al. (2002) we obtain, respectively, $`(mM)_{0,I}`$=18.49$`\pm `$0.06 and $`(mM)_{0,K}`$=18.46$`\pm `$0.03 (the errors include both the error on the reddening and the error on the local RC magnitudes). For the SMC, again using the observed values tabulated by PGU03, we obtain $`(mM)_{0,I}`$=18.90$`\pm `$0.06 and $`(mM)_{0,K}`$=18.88$`\pm `$0.03. The difference in the RC distance moduli, $`\mathrm{\Delta }M_{LMCSMC}`$, is equal to 0.41$`\pm `$0.08 in $`I`$ and 0.42$`\pm `$0.04 in $`K`$, which compares well with the value of 0.50$`\pm `$0.08 obtained from the comparison of the Tip of the Red Giant Branch magnitudes found by PGU03.

Significantly, the LMC distance moduli we find above, derived individually from the $`I`$\- and $`K`$-band data, are fully in agreement with the LMC distance modulus determined from the multiwavelength method of Alves et al. (2002). This method imposes the constraint that the derived distance modulus must be the same in all the passbands used – solving simultaneously in several passbands then yields both distance modulus and reddening. Using $`V`$, $`I`$ and $`K`$ data, Alves et al. (2002) find $`(mM)_{0,LMC}=18.493\pm 0.033_\mathrm{r}\pm 0.03_\mathrm{s}`$. Most recently, Salaris et al. (2003) have applied the same method to the field surrounding the LMC cluster NGC 1866 and found $`(mM)_{0,LMC}=18.53\pm 0.07`$ from $`HST`$ $`V`$ and $`I`$ data.

## Acknowledgments

This research is based on observations made at the South African Astronomical Observatory.

We have made extensive use of the WEBDA database (http://obswww.unige.ch/webda). We thank the anonymous referee for a prompt and careful reading of the manuscript. SMP acknowledges financial support from PPARC.


text:
# Estimating good discrete partitions from observed data: symbolic false nearest neighbors

(April 2, 2002)

## Abstract

A symbolic analysis of observed time series data requires making a discrete partition of a continuous state space containing observations of the dynamics. A particular kind of partition, called “generating”, preserves all dynamical information of a deterministic map in the symbolic representation, but such partitions are not obvious beyond one dimension, and existing methods to find them require significant knowledge of the dynamical evolution operator or the spectrum of unstable periodic orbits. We introduce a statistic and algorithm to refine empirical partitions for symbolic state reconstruction. This method optimizes an essential property of a generating partition: avoiding topological degeneracies. It requires only the observed time series and is sensible even in the presence of noise when no truly generating partition is possible. Because of its resemblance to a geometrical statistic frequently used for reconstructing valid time-delay embeddings, we call the algorithm “symbolic false nearest neighbors”.

Why might one want to represent observed time series of dynamical systems as sequences of low-precision discrete symbols? In this representation, there are often interesting techniques—often (but not exclusively) derived from information theory and its associated technology—which may illuminate data in novel ways symbolictimeseriesanalysis . The initial step for all these methods requires making a partition: a coloring of the state space reconstruction , $`𝐱R^d`$, into non-overlapping regions and associated symbols so that any $`𝐱`$ is assigned a unique symbol $`s`$ in a discrete alphabet. The symbol may be represented as an integer in the set $`0,1,\mathrm{}A1`$. A partition $`𝒫`$ defines a discretization of the observed sequence $`𝐱_i,i=1\mathrm{}N`$ into a symbolic sequence, $`s_i,i=1\mathrm{}N`$.

What partitions are “good”? Which discretizations retain the full structure of the original dynamics in the $`𝐱`$ space in the sequence of symbols? Unfortunately, the situation is unlike the remarkable time-delay embedding method for continuous dynamics: simple partitions are not generically satisfactory. The mathematics of symbolic dynamics specifies what we want: a “generating partition” (GP), where symbolic orbits uniquely identify one continuous space orbit, and thus the symbolic dynamics is fully equivalent to the continuous space dynamics.

Unfortunately there is no satisfactory mathematical theory about how to find a GP as a general procedure (except for one dimensional dynamics ($`d=1`$), where partitioning at the critical points works). Are the ad-hoc partitions often used still satisfactory? Unfortunately they are often not so. Bollt et al bollt  examined the degradation in the symbolic dynamics which results from the frequently used “histogram partition”, as opposed to a GP. A less optimal partition will induce improper projections or degeneracies, where a given symbolic segment may correspond to more than one topologically distinct state space orbit. This resulted in finding the wrong topological entropy. Chaotic communication with symbolic targeting works most satisfactorily knowing a GP, because then the transmitted symbolic message may be directly mapped into a desired orbit in the attractor. (see, e.g. lai2dsaddle )

Since a partition is a critical first step for any symbolic data analysis, a poor partition yield poor results, a method to approximate good partitions from observed data alone is urgently needed. With apparently no existing satisfactory solutions, this is the problem we attack. Davidchack et al davidchack  recently presented a partitioning method which works by successively coloring unstable periodic orbits (UPOs) to ensure unique codings (all UPOs have unique codes under a GP). The necessary high-order UPOs are very difficult to obtain from observed data alone, unfortunately.

The Kolmogorov-Sinai entropy rate $`h_{KS}`$ of the dynamics can be found from the supremum, over all increasingly fine partitions, of Shannon’s entropy rate evaluated on the information source implied by the discretization. More strikingly, a GP also achieves this supremum with a finite, and, one hopes, small alphabet. This suggests a naive strategy whereby one maximizes a statistical estimator of the entropy rate evaluated on the sequence induced by candidate partitions ebeling . This apparently attractive idea is flawed as demonstrated by the following counterexample. Consider a partition of the state space with a fine box-size $`ϵ`$ where each region is randomly assigned a symbol from the alphabet. For sufficiently small $`ϵ`$, the symbol sequence of any finite time series will appear indistinguishable from a memoryless and structureless information source with the maximum possible entropy rate $`h=\mathrm{log}_2A`$, since each observed datum could have encountered a different partition element with a new random symbol.. With the highest entropy possible this partition would be selected over competitors but is clearly useless for data analysis as the resulting symbolic stream says nothing about the original time series or its particular dynamics. Even if one believes this pathology to be irrelevant to coarser encodings, there are practical problems with the maximum-estimated-entropy idea. First, estimation of $`h`$ is not trivial to do well; second, when there is observational noise (inevitable with data acquisition equipment) larger alphabets will appear to give significantly higher entropies even if they are not actually much better at encoding the dynamics. As the true entropy rate of the system is not already known (and often a key quantity one wants to estimate given a good partition), there is no absolute statistical target which confirms whether the proposed partition is at all close or far from the ideal. In practice, selecting partitions with entropy does not seem to work well in general.

We assert our practical criterion for a good partition: short sequences of consecutive symbols ought to localize the corresponding continuous state space point as well as possible. A good coding ought to maintain the benefits of a low-precision symbolic representation with minimum distortion of the original state space dynamics. Our central idea is to form a particular geometrical embedding of the symbolic sequence under the candidate partition and evaluate, and minimize, a statistic which quantifies the apparent errors in localizing state space points.

We embed the symbol sequence into the unit square henonsymbols :

$$𝐲_i=(\underset{k=1}{\overset{k_{\text{max}}}{}}s_{i(k1)}/A^k,\underset{k=1}{\overset{k_{\text{max}}}{}}s_{i+k}/A^k).$$
(1)
($`k_{\text{max}}`$ is chosen such that $`A^{k_{\text{max}}}`$ is as small as the computational precision.) For a binary alphabet $`(A=2)`$, the first coordinate of $`𝐲_i`$ is the binary fraction whose digits start at $`s_i`$ and go backwards in time, the second is with the sequence going forward from $`s_{i+1}`$. Intuitively, the distribution on $`𝐲`$ is like a $`𝒫`$-dependent symbolic version of the invariant measure.

Given $`𝐱_i`$ and a partition $`𝒫`$, the symbolic embedding (1) yields a parallel series $`𝐲_i`$, defining points on some map $`𝐲=\varphi _𝒫(𝐱)`$. We want $`\varphi _𝒫`$ to be injective, i.e. $`\varphi _𝒫(𝐱)=\varphi _𝒫(𝐱^{})`$ implies $`𝐱=𝐱^{}`$. With finite data, we desire that if $`\varphi _𝒫(𝐱)\varphi _𝒫(𝐱^{})`$ is small, so is $`𝐱𝐱^{}`$. By construction, sufficiently near points in $`𝐱`$ have close symbolic sequences in their most significant digits. In a good partition, additionally, nearby points in $`𝐲`$ remain close when mapped back into the $`𝐱`$-space. By contrast, bad partitions induce topological degeneracies where similar symbolic words map back to globally distant regions of state space. As shown in bollt , this phenomenon confounds proper analysis of the observed symbolic dynamics.

We need to quantify how well any candidate partition achieves our ideal. We find the nearest neighbor, in Euclidean distance, to each point $`𝐲_i`$. Conventional $`k`$-d tree algorithms kdtree  efficiently provide the index of the nearest neighbor to any point in a data set: $`𝒩[i]=\text{arg}\mathrm{min}_{ki}𝐲_k𝐲_i`$. Knowing symbolic neighbors, we find distances of those same points back in $`𝐱`$-space, $`D_i=𝐱_{𝒩[i]}𝐱_i`$. We normalize the set of $`D_i`$ by a monotonic transformation: given any $`D`$, find its rank $`R[0,1]`$ in the cumulative distribution of random two-point distances $`𝐱_\alpha 𝐱_\beta `$. Large $`R`$ means that localizing well in symbol space did not localize well in the original state space.

Better partitions give a smaller proportion of symbolic false nearest neighbors, that fraction of $`R_i`$ which are greater than some threshold $`\eta `$, denoted $`J_{\text{sfnn}}`$. This resembles the false neighbors statistic for time-delay embeddings falseneighbors : both count large-deviation “mistakes” in a related space which result from topological mis-embedding in the tested space. Appropriate values for $`\eta `$ which defining a “large” deviation are $`\eta 0.010.3`$, depending on the noise in $`𝐱_i`$. An alternative to $`J_{\text{sfnn}}`$ is $`K_{\text{sfnn}}`$, defined as the arithmetic average of the largest $`\gamma `$ percentile of the set of $`R_i`$. Using $`J_{\text{sfnn}}`$, $`\eta `$ may need tuning depending on the noise scale and dynamical system, but the effect of changing $`\gamma `$ is lower. On the downside, $`K_{\text{sfnn}}`$ does not necessarily converge to near zero for the optimal partition. We typically find good results with $`\gamma 0.010.05`$.

For concrete numerical calculations, we need to parameterize partitions with a relatively small number of free parameters. Inspired by davidchack , we define partitions with respect to a set of radial-basis “influence” functions of the form $`f_k(𝐱)=\alpha _k/𝐱𝐳_k^2`$, the set of $`\alpha `$ and $`𝐳`$ being the free variables. For any particular $`𝐱`$, one $`f_l(𝐱)`$ will generically result in the largest value versus other $`f_k(𝐱),kl`$, and then $`𝐱`$ is assigned to that symbol which was pre-assigned to influence function $`f_l`$. The $`𝐳`$ parameters are initialized to random examples from the $`𝐱_i`$ and $`\alpha `$ to independent random variates $`[0,1)`$, and $`n_f`$ functions assigned to each of the $`A`$ symbols. In davidchack  the $`𝐳_k`$ were fixed on the UPOs and their symbols varied; here, the centers and coefficients vary but their symbols are fixed. We minimize $`J_{\text{sfnn}}`$ or $`K_{\text{sfnn}}`$ over the $`An_f(d+1)`$ free parameters using “differential evolution” diffevol , a genetic algorithm suitable for continuous parameter spaces.

Figure 1 shows the final $`𝒫`$ on 2000 data points from the Ikeda mapikedamap . It shows the best result (lowest $`J_{\text{sfnn}}`$) out of six restarts changing only the random seed governing the initial conditions; the results were not much worse on the other runs, however. The result is very close to the partition knowing the dynamics.

In a stationary information source, the number of distinct length-$`p`$ codewords will scale, for asymptotic $`p`$, as $`N_pe^{h_Tp}`$ where $`h_T`$ is the topological entropy, a dynamical invariant. We validate $`𝒫`$ with an estimate of the deficiency between $`h_T`$ implied by $`𝒫`$ and the correct $`h_T`$: $`\delta h_T=p_{\text{max}}^{}{}_{}{}^{1}_{p=1}^{p_{\text{max}}}p^1\mathrm{log}\left(N_p/\stackrel{~}{N}_p\right)`$. $`N_p`$ is the number of distinct period-$`p`$ UPOs (which were computed knowing the equations of motion), $`\stackrel{~}{N}_p`$ the number of such UPOs with unique $`p`$-symbol codes in some $`𝒫`$. A GP gives $`\delta h_T=0`$, and $`\delta h_T0`$ for better (less UPO-degenerate) partitions. Figure 2 shows $`\delta h_T`$ on each new best partition found during the optimization. The optimization target, $`J_{\text{sfnn}}`$, decreases strictly monotonically by construction; though $`\delta h_T`$ does not decrease quite monotonically, the trend toward very small values is clear. This gives evidence that minimizing $`J_{\text{sfnn}}`$ also refines approximations to GPs.

Figures 35 demonstrate applications of the algorithm. Fig. 3 shows the effect of noise on a system where the GP is analytically known. The Lozi map (see analysis in henonsymbols ) is similar to the Henon map but replaces the quadratic nonlinearity with a piecewise linear one. We find a partition which is close to the noise-free GP even when the data have been contaminated by significant amounts of additive noise. Though complete localization to a single point is not possible here, minimizing large divergences is still a desirable criterion. Figures 4 to 5 show estimated partitions on experimental data sets where no analytical form of the equations (much less partitions) are known. On account of noise, dynamical or observational, a certain amount of divergence $`D`$ for symbolic nearest neighbors is inevitable. Still, minimizing large deviations is a reasonable goal even for noisy data. There are very few rank distances with $`R0.2`$ or $`0.3`$ compared to a basic histogram partition.

It must be kept in mind that GPs are not necessarily unique for any attractor. Distinctly different partitions may be found, all of which are reasonably satisfactory. (At a minimum any iterate of a GP is also a GP). There are many coexisting solutions (roughly like finding low-energy states of a spin glass), which is why the optimization problem is hard, requiring a global search method. We conjecture this is one reason why understanding the structure of generating partitions has been so difficult for mathematicians.

We also point out that it is possible to partition one space $`𝐱`$, but quantify distances $`D_i=𝐳_{𝒩[i]}𝐳_i`$ in some different space as long as there is some relation between each $`𝐱_i`$ and $`𝐳_i`$. For example, one may be interested in a simple symbolic control scheme, say where $`𝐳_i=𝐱_{i+T}`$ (find the best partition of observables now that best predicts some future), or perhaps when two different variables are measured simultaneously and one wants to cross-predict. In the first case, in principle, a GP should be optimal for predicting the future as well as the present but the inevitable issues of finite data and noise may make the best empirical partition different for the two cases. In the second case, generating partitions are irrelevant entirely.