text:
# Quantum fluctuations and correlations of spatial scalar or multimode vector solitons in Kerr media

## 1 Introduction

The analysis of spatial aspects of quantum fluctuations is now a well established field of quantum optics (see for example  for a review). More specifically, some recent studies  have considered quantum effects in spatial solitons. Ref.  shows that vacuum-induced jitter remains negligible versus the width of a soliton. In Ref. , production of subpoissonian light is demonstrated by placing either a stop on the centre of the beam in the near field or a low pass filter in the far-field. These two theoretical papers use linearisation of quantum fluctuations. Ref.  gives a general method to numerically determine the transverse distribution of linearised quantum fluctuations by using Green functions. Results show an important global squeezing for $`\chi ^2`$ as well as for $`\chi ^3`$ solitons. In Kerr media however, this squeezing does not occur on the amplitude quadrature, implying homodyne detection. Ref. demonstrates also local squeezing and anti-correlations between different parts of the beam. The Green’s function method has also been applied to the assessment of spatial aspects of quantum fluctuations of spontaneous down-conversion . In this paper, we apply the Green’s function method to spatial scalar solitons and to multimode vector solitons in Kerr media . The results for the scalar solitons appear somewhat different of that given in ref. , because of a flaw in the numerics used in this reference. We show in particular that, in the corrected results, the squeezing monotonically increases with the propagation length, like for a plane wave. We show also that results of ref.  are retrieved with the Green’s function method.

We study also the quantum properties of multimode vector solitons that result from incoherent coupling of two counter-rotating circularly polarised fields mutually trapped through cross-phase modulation, the self-focusing nature of a component being exactly compensated by the diffracting nature of the other. We find the quantum limits for the stability of these multimode vector solitons and show that the global squeezing (measured on the total intensity of the soliton) appears to be similar to that of scalar solitons. In addition, we demonstrate that there exists a strong anti-correlation between the circular components.

The paper is organised as follows. In Sec. 2, we recall the principle of the Green’s function method. Sec. 3 gives the results for scalar Kerr solitons and Sec. 4 is devoted to vector solitons.

## 2 Green’s function method

Let us briefly recall the main features of the Green’s function method  that allows us to determine the spatial distribution of quantum fluctuations in the regime of free propagation in a non-linear medium. We will present it in the simplest case of a scalar field and in the configuration of a single transverse coordinate x. The classical propagation equation for the complex electric field envelope $`U(x,z)`$ is in this case :

$$\frac{U}{z}=\frac{i}{2k}\frac{^2U}{x^2}+i\gamma |U|^2U$$
(1)
where $`\gamma `$ is the usual Kerr coefficient, z the main propagation direction, $`k`$ the longitudinal wavevector, and x the transverse position in the transverse plane. This equation has a solution that we will call $`\overline{U}(x,z)`$.

Let us write the quantum positive frequency operator describing the same field at the quantum level, $`\widehat{U}(x,z)`$, as :

$$\widehat{U}(x,z)=\overline{U}(x,z)+\delta \widehat{U}(x,z)$$
(2)
If the fluctuations remain small compared to the mean fields, which is true as long as the field contains a macroscopic number of photons, then the quantum mean field coincides with the field given by classical nonlinear optics, $`\overline{U}(x,z)`$, and the fluctuations $`\delta \widehat{U}`$ obey a simple propagation equation, obtained by linearizing the classical equation of propagation 2 around the mean field. This equation writes :

$$\frac{}{z}\delta \widehat{U}=\frac{i}{2k}\frac{^2}{x^2}\delta \widehat{U}+i\gamma \left(2|\overline{U}|^2\delta \widehat{U}+\overline{U}^2\delta \widehat{U}^+\right)$$
(3)
As this equation is linear, one can use Green’s function techniques to solve it. Its solution at a given $`z_{out}`$ is a linear combination of the input fluctuations operators in plane $`z_{in}`$:

$$\delta \widehat{U}(x,z_{out})=𝑑x^{}G(x,x^{})\delta \widehat{U}(x^{},z_{in})+𝑑x^{}H(x,x^{})\delta \widehat{U}^+(x^{},z_{in})$$
(4)
where $`G(x,x^{})`$ and $`H(x,x^{})`$ are the Green’s functions associated with the linear propagation equation 3. If the input field fluctuations are those of a vacuum field, or of a coherent field, they are such that :

$$<\delta \widehat{U}(x,z_{in})\delta \widehat{U}^+(x^{},z_{in})>=C\delta (xx^{})$$
(5)
where $`C`$ is a constant. Equations 4 and 5 allow us to calculate quantum noise variance or correlation at the output of our system on any quadrature component when one knows the two Green’s function $`G(x,x^{})`$ and $`H(x,x^{})`$ of the problem. These two functions can readily be numerically calculated as solutions of the propagation equation 3 when the initial condition is a delta function localised at a given transverse point. This method can be obviously extended to any kind of non-linear propagation equation.

We will focus our attention here on two kinds of simple spatially resolved measurements : the first one is the direct photodetection on a pixellised detector, which gives access to the local intensity fluctuations ; the second one is the spatially resolved balanced homodyne detection, using a local oscillator having a plane wavefront, coincident with the detector plane and with a global adjustable phase $`\theta `$, which allows us to measure the local fluctuations $`\delta \widehat{I}(\theta ,x)`$ of a given quadrature component. Reference  gives the mathematical expressions, in terms of the Green’s functions $`G(x,x^{})`$ and $`H(x,x^{})`$, of these quantities, from which one easily derives the best local squeezing $`min_\theta <(\delta \widehat{I}(\theta ,x))^2>`$ and the covariance between pixels :

$$C(x,x^{},\theta )=<\delta \widehat{I}(\theta ,x)\delta \widehat{I}(\theta ,x^{})>$$
(6)
## 3 Results : scalar solitons

Equation (1) describes the propagation of the electric field envelope of a monochromatic field of frequency $`\omega `$ in a planar waveguide made of a Kerr medium, i.e. having an intensity dependent index of refraction:

$$n=n_0+n_2|U|^2$$
(7)
with $`k=n_0\omega /c`$ and $`\gamma =n_2\omega /c`$.

It can be written in a universal form if one uses a scaling parameter $`\eta `$ and the scaled variables $`\zeta =\eta z`$, $`r=x\sqrt{2\eta k}`$, and its solution is the well known $`\zeta `$ invariant hyperbolic secant proportional to $`1/coshr`$ : the 1D scalar soliton, which is able to propagate without deformation in the nonlinear medium.

### 3.1 Squeezing on the total beam

As it is well known , the Kerr effect produces in the plane-wave case a significant amount of squeezing, which increases monotonously with the propagation distance. The squeezing is optimized for a given quadrature component (best squeezing) which is neither the amplitude nor the phase quadrature. However, the amplitude quadrature of the field remains at the shot noise level. Fig 1 gives the results of the Green’s function method for the scalar soliton, and shows that the results are quite similar to the plane wave case, in the case of a photodetector having an area much larger than the soliton spot. Neglecting diffraction, or using single mode propagation in an optical fiber  seems almost equivalent to compensating diffraction by self-focusing, for a given nonlinear phase $`\zeta `$. Hence a spatial soliton appears as a practical mean to obtain strong squeezing, with the restriction that it must be detected by homodyne techniques.

### 3.2 Local quantum fluctuations.

In order to measure the local quantum fluctuations, we use a homodyne detector, made of photodetectors of very small area and quantum efficiency 1, and using a local oscillator at frequency $`\omega `$, placed at the output of the nonlinear medium, which is able to monitor the quantum noise on any quadrature component at point $`(x,z_{out})`$ and we vary the local oscillator phase to reach the minimum noise level : we obtain by this way a quantity that we call ”best squeezing”. Figures (2a) and (2b) give, for two normalized propagation distances ($`\zeta =0.3`$ and $`\zeta =3`$), such a quantity as a function of x, together with the phase angle of the local oscillator enabling us to reach the best squeezed quadrature. One observes that the central pixel is the most squeezed. However, the noise level appears to be above the shot noise for the longer propagation length.

Indeed, figure 3, which displays the variation of the best squeezing on the central pixel as a function of propagation length, shows that the fluctuations on this pixel decrease until an optimal propagation length ($`\zeta =0.3`$), then increase and overpass the shot noise for a sufficiently long propagation length. Hence, long propagation lengths produce at the same time a local excess noise and a squeezed total beam (see fig.1). These results are not contradictory because, as we will soon see, there is a gradual build-up of anti-correlations between the different transverse parts of the soliton because of diffraction.

Figure 4 shows that the covariance between pixels on the best squeezed quadrature (calculated with eq.6 and where best squeezing is defined for the total beam) increases and spreads out when passing from $`\zeta =0.3`$ to $`\zeta =3`$. Indeed, only pixels close to each other are anti-correlated for $`\zeta =0.3`$ while the correlation between adjacent pixels becomes positive for $`\zeta =3`$. For this latter distance, squeezing is due to the anti-correlation between the left and the right side of the soliton.

Figure 5 shows the best squeezing for a photodetector of variable size, centered on the beam. This photodetector can be seen as an iris that allows the detection of only the central part of the beam, with a minimum size corresponding to the central pixel and a maximum size corresponding to the detection of the entire beam. As expected, the noise on a very small area is close to the shot noise. The squeezing on a small central area is bigger at small distance, while the squeezing for the entire beam is better at great distance. For $`\zeta =0.3`$ the squeezing is maximum for a finite size of the photodetector (transmission coefficient of 0.8, corresponding to a normalized radius of 1), while detecting the entire beam ensures the maximum squeezing for $`\zeta =3`$.

In many circumstances, a spatial soliton behaves as a single mode object. The present analysis allows us to check such a single mode character at the quantum level : let us assume that the system is in a single mode quantum state. This means that it is described by the state vector $`|\mathrm{\Psi }>|0>\mathrm{}|0>\mathrm{}`$, $`|\mathrm{\Psi }>`$ being a quantum state of the mode having the exact spatial variation of the soliton field, and all the other modes being in the vacuum state. It has been shown in  that if a light beam is described by such a vector, then the noise recorded on a large photodetector with an iris of variable size in front of it has a linear variation with the iris transmission. We see in figure 5 that it is not strictly speaking the case : as far as its quantum noise distribution is concerned, a spatial soliton is not a single mode object. Let us also notice that the departure from the linear variation is on the order of $`10\%`$, so that the single mode approach of the soliton remains a good approximation.

### 3.3 Intensity squeezing by spatial filtering

In all results presented in the preceding paragraph, squeezing can be measured only by homodyne detection, because no squeezing is present for the amplitude quadrature, which is the easiest to measure, as it requires only direct detection, and no interferometer. Actually, Mecozzi and Kumar  have shown that a simple stop on the central part of the soliton beam ensures intensity squeezing on the remaining of the light. They have proposed the following intuitive explanation : if the intensity fluctuation is positive for the entire beam, the soliton becomes narrower (its width is inversely proportional to its power) and the stop will produce larger losses. In the Fourier plane, a square diaphragm (low-pass filter) will produce an even better squeezing. Figure 6 shows the intensity squeezing we have obtained with the Green’s function method by considering an aperture placed in the Fourier plane. These results are in full agreement with the fig. 1 of Ref. 3, that uses another linearisation method. The residual differences come probably from a small difference in the aperture width. To conclude this paragraph, we find, as Mecozzi and Kumar, that subpoissonian light can be very easily produced from a spatial soliton by placing a simple aperture in the Fourier plane.

## 4 Quantum limits of stability and squeezing on vector solitons

### 4.1 Multimode vector solitons and symmetry-breaking instability

In the last few years, a great number of multi-component solitons have been demonstrated . Among these, the simplest example is the bimodal vector soliton which is depicted in fig. 7a and consists of two opposite circular components that trap each other in a non birefringent Kerr medium. Their complex envelopes U and V obey the following coupled nonlinear Schrödinger equations :

$$\{\begin{array}{cc}\frac{U}{Z}=i\frac{1}{2k}\frac{^2U}{x^2}+i\gamma (|U|^2U+7|V|^2U)\hfill & \\ \frac{V}{Z}=i\frac{1}{2k}\frac{^2V}{x^2}+i\gamma (|V|^2V+7|U|^2V)\hfill & \end{array}$$
(8)
Where $`\gamma `$ is the usual Kerr coefficient, and 7 represents the strength of the cross-phase modulation in an isotropic liquid, like $`CS_2`$ . As can be seen on fig.7b, this soliton bound-state propagates in equilibrium due to incoherent coupling between both circular components. If alone in the Kerr medium, the symmetrical field U tends to self-focus. However, the antisymmetrical field V on the counter-rotating polarization tends to diffract and it can be shown that a proper choice of the intensities (in dimensionless units on Figs 7) ensures an exact balance between diffraction and self-focusing. However, this equilibrium is unstable if the cross-phase modulation coefficient is greater than the self-phase modulation coefficient. Figure 8a displays the solution of the classical equations of propagation when a very small spatial initial asymmetry is introduced in one of the two fields. It shows that the multimode soliton is unstable and eventually collapses, the intensity evenly distributed between the two cores of the induced waveguide going abruptly from one core to the other. The instability occurs after a propagation distance that decreases as the initial asymmetry increases, as it can be clearly seen in fig. 8b. A more complete description of this phenomenon is given in Ref. , while the first experiment demonstrating symmetry-breaking instability is reported in Ref. . In this experiment, performed with pulses issued from Nd-Yag at a 10 Hz repetition rate, symmetry-breaking induced by the noise at the input occurred for about $`60\%`$ of the laser shots, with a left-right repartition equally and randomly distributed.

### 4.2 Stability of vector solitons : quantum limits

Propagation of a multimode vector soliton in a medium with strong cross phase-modulation is unstable, as shown theoretically by Kockaert et al. , using a linear stability analysis. Because of quantum noise, symmetry-breaking will occur even for an ideal symmetrical input beam. The simplest way to quantitatively assess the phenomenon is propagating a vector soliton, using Eq.(8) and the usual split-step Fourier method, after the addition on the input fields of a Gaussian white noise with a variance of half a photon per pixel on each polarization. The validity of such type of stochastic simulation is well known in the framework of the Wigner formalism . Results show symmetry-breaking with a random direction after 60 to 70 diffraction lengths.

A linear stability analysis  shows that exponential amplification occurs only for the part of the input noise obtained by projection on the eigenvectors associated to eigenvalues with a negative imaginary part. Such an eigenvector corresponds physically to a transverse distribution of field for each circular polarization. The part of an unstable eigenvector that is circularly polarized as the U symmetrical field is antisymmetrical, while the part of the eigenvector that is circularly polarized as the V antisymmetrical field is symmetrical. Moreover, exponential amplification corresponding to the eigenvalue having the most negative imaginary part is predominant after a short propagation distance. Hence, instability can be described as the amplification of a single unstable eigenvector. In the results of figure 9, the field corresponding to this eigenvector with an energy of half a photon per polarization has been added at the input to the ideal vector soliton. Figure 9 shows in solid line the evolution of the norm of the difference between the input field and the propagated field, compared to a pure exponential amplification with the corresponding imaginary eigenvalue (dashed line). The difference evolves almost exponentially, with however a periodic modulation due to the real part of the eigenvalue. This phase modulation leads to a dependence between the direction of symmetry-breaking and the initial amplitude of the perturbation. Symmetry-breaking appears clearly on figure 9 when the initial small perturbation has grown until a macroscopic level, for a propagation distance of about 70 diffraction lengths.

To conclude this paragraph, the usual classical simulation of quantum noise allows us to predict a maximum propagation length for the multimode vector soliton. This length is less than two meters, for the diffraction length characterizing the experiment of ref., but also 70 times the typical length leading to symmetry-breaking in this experiment. Actually, experimental evidence of quantum limits for symmetry-breaking appears as a rather difficult challenge, because unavoidable spatial classical noise dominates in realistic input fields.

### 4.3 Total beam squeezing and correlation between polarisations

In this paragraph, we give the results of the Green’s function method in the case of vector solitons. Because these solitons are not analytically integrable, the impulse response due to a Dirac-like field modification of the perfect vector soliton is more complicated to determine, and calculated as follows. The perturbation consists in a unity single-pixel step, multiplied by a small coefficient before to be added to the field in order to ensure a near-perfect linearity of propagation equations for the perturbation. Both the unmodified and the modified field are first numerically propagated using eq.8, then subtracted each other. Finally, the subtraction result is divided by the initial multiplication coefficient, in order to retrieve the output corresponding to the input unity single-pixel step. An important point is that the phase reference corresponding to the amplitude quadrature is given for each pixel by the unmodified field. Indeed, in practice the vector soliton could experience a very slight breathing, that results in a weak phase curvature. Therefore, taking the mean soliton phase as reference leads to strong numerical inaccuracies for strong squeezing.

Fig 10 shows the global squeezing coefficient of the vector soliton on the best quadrature. This figure is very similar to that obtained for a scalar soliton (figure 1) : strong squeezing occurs, although somewhat smaller than in the scalar case for the same propagation distance. However, this squeezing is due to a strong anti-correlation between circular components, as shown in figure 10 b. It can be verified on figure 11, which displays the best squeezing of each circular polarisation component as a function of propagation, that the best squeezing on a single circular polarization is much weaker. We can conclude that a vector soliton as a whole exhibits quantum properties that are similar to a scalar soliton, while the field on each polarisation does not constitute itself a soliton.

### 4.4 Local quantum fluctuations

Figure 12 shows the squeezing coefficient on one pixel versus its transverse position, for the propagation length where this local squeezing is the most intense. The degree of squeezing is similar to the scalar case (see figure 2), with however two peaks of squeezing located on the intensity peaks of the multimode vector soliton. As in the scalar case, local squeezing is maximum for a relatively small propagation length (fig.12 b). Figure 13 shows the best squeezing for a photodetector of variable size, centered on one of the intensity peaks. In contrast to figure (5) the variation is far from being linear with the transmission, meaning that the system is not single mode, as expected. One observes that for long propagation distances, the global squeezing reaches a maximum for a photodetector size corresponding roughly to the peak width. This global squeezing implies the existence of spatial anti-correlations. Indeed, the local squeezing on one pixel disappears for such long distances (see figures 12 and 13).

The correlations between fluctuations at different pixels, and possibly on different polarizations, provide a great deal of information about the distribution of fluctuations in the soliton. Figure 14 shows these spatial correlations, in terms of the covariance functions $`C(x,x^{},\theta )`$ of the two polarizations, and the covariance function between the two modes, on the best squeezed quadrature and for 2 propagation lengths. One sees that for $`\zeta =0.6Ld`$, anti-correlations appear between pixels when one measures different circular polarizations on the pixels. It means that the effect of cross-phase modulation is already noticeable, while effects of diffraction and self-phase modulation are still weak. For $`\zeta =2Ld`$, anticorrelations appear also between pixels on the same polarization, as for the scalar soliton.

### 4.5 Fluctuations on linear polarisations

Figure 15 shows that squeezing on linear polarisations disappears after a relatively short propagation distance, even on the best quadrature. Moreover, fluctuations grow exponentially. On the other hand, fluctuations on both linear polarisations are perfectly anti-correlated after some distance, as it was expected because of the strong squeezing of the global beam. We have verified that the results on this global beam (section 4.3) can be retrieved also from a description of the field by projection on linear polarisations.

## 5 Conclusion

Let us summarise here the main results about the spatial distribution of quantum fluctuations in 1D Kerr solitons that we have obtained in this paper : as far as measurements on the whole soliton are considered, we have shown that scalar Kerr solitons, as well as vector Kerr solitons, exhibit squeezing properties similar to that of a plane wave propagating in the same medium and for the same type of propagation distance. Because of diffraction, as the soliton propagates, strong anti-correlations develop between symmetrical points inside the scalar soliton spot, and between polarisation components for the vector soliton.

## 6 Acknowledgments

This work has been supported by the European Union in the frame of the QUANTIM network (contract IST 2000-26019).

## 7 References


format: text

text:
# Efficiency and reliability of epidemic data dissemination in complex networks

## Abstract

We study the dynamics of epidemic spreading processes aimed at spontaneous dissemination of information updates in populations with complex connectivity patterns. The influence of the topological structure of the network in these processes is studied by analyzing the behavior of several global parameters such as reliability, efficiency and load. Large scale numerical simulations of update spreading processes show that while networks with homogeneous connectivity patterns permit a higher reliability, scale-free topologies allow for a better efficiency.

Modern society increasingly relies on large scale computer and communication networks, such as the Internet. A major challenge in these networks is the developing of reliable algorithms for the dissemination of information from a given source to thousands, or even millions, of users such as for news and stock exchange updates, mass file transfers and Internet broadcasts dave ; deering . In epidemic-inspired communication this is achieved by exploring a mechanism analogous to the spreading of infectious diseases in populations vogels ; ep3 . Indeed, epidemic data dissemination in computer and communication networks show interesting parallels with both disease propagation in populations and the spread of rumor in social networks. The information spreads like a benign epidemic through local interaction between nodes which forward the message they receive to a random selection of their peers in the network, until the whole system becomes “infected” with information. The great advantages of epidemic-style communication is that dissemination proceeds on a local basis, without any co-ordination from a central organizing body vogels ; ep3 . These protocols are also highly resilient to sudden failure of communication links and nodes.

A relevant result in the mathematical theory of epidemics is that the spreading of infection in a population is strongly affected by the patterns of connectivity in the underlying contact networks. In particular in scale-free topologies, characterized by degree distributions with power law behavior $`P(k)k^\gamma `$ bar99 ; fal ; romu01 , the statistical relevance of hubs makes the network highly permeable to attacks newman00 ; havlin00 ; bar00  and the spreading of infections pv01a  and highlights the need of special immunization strategies. This result suggests that the topology of the underlying computer and communication network might heavily affect the performance of epidemic-style data dissemination protocols. Surprisingly, however, the impact of network topology on such protocols has not been thoroughly explored, although the results could have an important technological value. Indeed, these protocols potentially find a large spectrum of application such as mobile communication networks and, more recently, resource discovery in the so-called peer-to-peer systems built on top of the Internet p2p , and finally in grid computing grid .

In this paper we define a simple epidemic data dissemination model and perform a detailed numerical study of the dynamics of the information propagation in networks with diverse topological properties. Our basic model is a slightly modified version of the Daley and Kendall (DK) model dk64 ; zanette ; liu  and it can be considered as the simplest epidemic algorithm for the updating of distributed databases ep1 ; ep3 . We study the main relevant features of the model such as the reliability of the dissemination process and the amount of traffic generated by the dynamics. Specifically, we study two different prototypical networks: a random homogeneous network note1  and a scale-free network. The results obtained point out that in the homogeneous topology the epidemic process provides a more reliable updating of the network. The scale-free topology, on the other hand, allows the algorithm to perform more efficiently in terms of the generated traffic. Finally, we compare the present model with a deterministic broadcast process and find that in a wide range of the model parameters the epidemic algorithm is more efficient.

The model we shall consider is defined in the following way. Each of the $`N`$ elements of the network can be in three possible states. We call a node holding an update and willing to transmit it a spreader. Nodes that are unaware of the update will be called ignorants while those that already know it but are not willing to spread the update anymore are called stiflers. We denote the density of ignorants, spreaders and stiflers at time $`t`$ as $`\psi (t)`$, $`\varphi (t)`$ and $`s(t)`$, respectively such that for all $`t`$, $`\psi (t)+\varphi (t)+s(t)=1`$. The spreading process takes place along the links between spreaders and ignorants. Each time step spreaders contact one (or more) neighboring node. When the spreader contacts an ignorant, the last one turns into a new spreader at a rate $`\lambda `$. On the other hand, the spreader becomes a stifler with rate $`1/\alpha `$ if a contact with another spreader or a stifler takes place note2 . The parameter $`\alpha `$ can be considered as the average number of contacts with spreader/stifler nodes before the spreader turns itself into a stifler. This dynamics mimics the attempt of diffusing an update or rumor by nodes which have been recently updated. At the same time, if a node attempts too many times to communicate the update to nodes which have already received it, it stops the process turning itself into a stifler. In other words, the node realizes that the update has lost its novelty and becomes uninterested in diffusing it. The present dynamics thus introduces a trade-off in maximizing the number of updated nodes and minimizing the number of contacts attempted. Obviously, the efficiency of the spreading process will depend on the rate at which individuals loose interest in further spreading of the rumor and the topology of the underlying network.

At the mean-field level and using the homogeneous assumption for the network connectivity pattern, the time evolution of ignorants, spreaders and stiflers is described by the simple set of equations

$$_t\psi (t)=\lambda \psi (t)\varphi (t),$$
(1)
$$_t\varphi (t)=+\lambda \psi (t)\varphi (t)\frac{1}{\alpha }\varphi (t)(\varphi (t)+s(t)),$$
(2)
where $`s(t)`$ is obtained by the normalization condition $`s(t)=1\psi (t)\varphi (t)`$ dk64 ; zanette . The dynamics asymptotically evolves to the state $`\varphi (\mathrm{})=0`$ in which the system is frozen. Noticeably, in random homogeneous networks, the density $`s(\mathrm{})`$ of elements which are aware of the update is always a finite fraction of the whole population dk64 ; zanette . The homogeneous assumption is, however, not valid anymore in the case of heterogeneous scale-free networks where it is known that spreading processes may show very different properties pv01a . In particular, an explicit dependence on the nodes’ degree $`k`$ must be included in the rate equations. While a general analytical solution cannot be obtained in this case, numerical studies on scale-free networks can be used to evaluate the reliability and efficiency of this process in more complex topologieszanette ; liu .

In the present investigation we used two specific network models. First we consider the Barabási-Albert (BA) network bar99 . In this model, starting from a set of $`m_0`$ nodes, one preferentially attaches each time step a newly introduced node to $`m`$ older nodes. The procedure is repeated many times and a network with a power law degree distribution $`P(k)k^\gamma `$ with $`\gamma =3`$ and average connectivity $`k=2m`$ builds up. This network is a clear example of highly heterogenous network, in that the degree distribution has unbounded fluctuations. As a reference of homogeneous networks we considered the Watts-Strogatz (WS) network ws98  in the case of complete random rewiring. In this case, one starts from a ring with $`N`$ nodes, each of them connected symmetrically to $`2K`$ neighbors. With probability $`p`$ each link connected to a clockwise neighbor is rewired to a randomly chosen node; otherwise it is preserved. After enough iterations a random network with an exponential connectivity decay for large $`k`$ and $`k=2K`$ is generated. Henceforth, we will use $`m_0=m=3`$ for the BA network and $`p=1`$ and $`K=3`$ for the WS model giving $`k=6`$ for both networks.

We have performed large scale numerical simulations by applying repeatedly the rules stated above on BA and WS networks. Initially, $`\psi (0)=\frac{N1}{N}`$, $`\varphi (0)=\frac{1}{N}`$, and $`s(0)=0`$, i.e., we start from a single spreader who is willing to spread the update through the network. At every time step, each of the $`\varphi N`$ spreaders contacts all its neighbors in a random sequence, unless during a contact it turns into a stifler. In this case it immediately stops contacting further nodes. This account for the larger transmission capabilities of high degree nodes that can reach a larger number of neighbors as specified by the heterogeneous network topology. The dynamical rules of the model are applied in parallel. The sizes of the networks used in the simulations carried out range from $`N=10^3`$ nodes to $`N=10^5`$ nodes and all numerical results have been obtained by averaging at least over 10 different networks and $`10^3`$ iterations. The parameter $`\lambda `$ may be varying as for the case of communication networks where it is known that the rate of packet loss is not always zero. Nevertheless, without loss of generality, $`\lambda =1`$ since it just fixes the time scale by rescaling opportunely the rate $`1/\alpha `$ in the Eq. 2. On the other hand, we vary the rate at which spreaders decide not to communicate the update any more from $`\alpha =1`$ to $`\alpha =10`$ and monitor several quantities of interest.

In order to characterize the propagation process we first focus on the reliability $`R`$ of the rumor propagation defined as the final density $`s(\mathrm{})`$ of nodes that have got the update when the process dies out. For obvious practical purposes, any algorithm or process that emulates an efficient spreading of a given message or data packet will try to raise as much as possible this magnitude.

In Table 1 we report the reliability of the spreading process for the BA network and the random graph generated with the WS algorithm with $`p=1`$ for several values of the parameter $`\alpha `$. As noticed previously liu , in the WS network the number of stiflers at the end of the process is already high even for $`\alpha =1`$, and the BA network appears less reliable. In general, it results that homogeneous networks allows a larger reliability $`R`$ to epidemic updating processes. This is not straightforward and one may think that the existence of hubs in scale-free networks helps propagate the rumor. However, a closer look at the spreading dynamics tells us that the presence of hubs introduces conflicting effects in the dynamics. While hubs may in principle reach a larger number of nodes, spreader-spreader and spreader-stifler interactions get favored on the long run. Indeed, it is very unlikely that a hub in the spreader state contacts all its ignorant neighbors before turning into a stifler. Once a few hubs are turned into stiflers many of the neighboring nodes could be isolated and never get the update. In this sense, homogeneous networks allow for a more capillary diffusion of the update, since all nodes contribute equally to the message-passing. This is opposite to what happens in the usual epidemic spreading model in heterogeneous networks models. These models lack of an infected recovery rate induced by neighbors already infected and fully exploit the advantage of hub’s large degree pv01a .

The previous discussion, however, refers only to average properties. In SF networks the connectivity distribution is highly heterogeneous and it is interesting to have a more detailed insight on the reliability of the process for different connectivity classes. It may be particularly relevant that a higher $`R`$ corresponds to the highly connected nodes, the hubs, which have a dominant role in the system. Figure 1 shows the behavior of the reliability $`R_k`$ measured as a function of the nodes connectivity $`k`$. This amounts to the relative density of nodes with connectivity $`k`$ that have received the update and it is measured as $`R_k=\frac{S_k}{N_k}`$, where $`S_k`$ and $`N_k`$ denote the total number of stiflers and the total number of nodes with degree $`k`$, respectively. The $``$ represents the average over many realizations. The results confirm that during the spreading dynamics it is very likely that highly connected nodes are reached by the update. In fact $`1R_k`$ decreases exponentially with $`k`$ and very high levels of reliability are obtained well before the natural cut-off of the network ($`k_{max}10^2`$) even for moderate values of $`\alpha `$. This is an interesting feature signaling that heterogeneous topologies can be considered reliable as far as the ”hubs” are concerned. We should note, however, that since the reliability is not a model independent quantity, this result might depend upon the specific details of the rumor algorithm.

In general, one does not only want to have high reliability levels, but also the lowest cost in terms of network load ep1 ; ep3 . This is generally achieved by imposing the minimum possible load to the network. Here, we define the load $`L`$ imposed to the network as the number of contacts established per node, i.e., how many messages on average each node sends to its neighbors in order to propagate the update. By using this quantity, an obvious definition of the global efficiency $`E`$ of the whole process is represented by the number of individuals who have got the update per unit of load, $`E=\frac{R}{L}`$. Its physical meaning is straightforward: the efficiency is equal to the fraction of “useful messages” (number of sites reached by the rumor) over total “load” imposed to the system. In Fig. 2, we report the behavior of the efficiency of the spreading process as a function of $`\alpha `$. In this case, the scale-free topology appears the most efficient for data dissemination. Indeed, the relative difference between SF and WS networks in global efficiency is larger than $`10\%`$ up to values of $`\alpha =5`$. This can be appreciated also by looking at the inset in Fig. 2. For both topologies the load on the network grows with $`\alpha `$, but the load imposed on SF networks is always smaller than on WS nets. Finally, it is interesting to compare the epidemic algorithm efficiency with those of the simplest broadcast strategy. This strategy essentially consists of a deterministic message-passing of each element to all its neighbors except the one from which the first update has been received. This way, a reliability $`R=1`$ is achieved since all nodes are surely contacted. In this case, the load is simply given by $`k1`$. In the case of networks with $`k=6`$ as the one used in the present study, the efficiency of the broadcast strategy is therefore $`E=0.2`$. It is interesting to note that in both the BA and the WS case, the epidemic algorithm achieve a better efficiency for a wide range of values of the parameter $`\alpha `$. This indicates that epidemic algorithms can provide attractive alternatives to broadcast solutions so far as the efficiency is concerned.

It is worth stressing that we are considering strategies in which the nodes do not have memory; i.e. they may try to resend a message to a node that has been contacted before. It is possible, however, to conceive different dynamics of the updating spreading strategies in which a trade-off between the memory introduced in the process and the optimization of reliability and efficiency is opportunely chosen. Other options rely on a careful tuning of the local message sending dynamics. For instance, we considered the case in which each spreader contacts only one node at each time step, reducing the effects of the hubs. This simple change allows for higher levels of efficiency, however, at the price of a much lower reliability. In general, thus it is possible to devise and tailor different processes that optimize one or more features of the update spreading on a given topology. We defer a more detailed study of this issue to future work.

In summary, we have studied the effect of the complex topological properties of many real networks in epidemic strategies for the communication of updates. The obtained results stimulate the seeking of heuristics and analytical methods to optimize epidemic algorithms taking into account the specific topology of the underlying network. These studies may have a large impact in technological and communication networks where the use of rumor-like algorithms might become a practice for data dissemination, reliable group communication or replicated database maintenance. Finally, they might provide a deeper understanding for social phenomena such as the spreading of new ideas in a population or the efficiency of marketing campaigns.

###### Acknowledgements.

We are grateful to A. Barrat and R. Pastor-Satorras for helpful comments and suggestions. Y. M. acknowledges financial support from the Secretaría de Estado de Educación y Universidades (Spain, SB2000-0357) and the hospitality of the Complexity Research Group of BT Exact, U.K. A. V. is supported by the EC - Fet Open project COSIN IST-2001-33555. This work has been partially supported by the Spanish DGICYT project BFM2002-01798.


text:
# Gaugino and meson condensates in 𝒩=1 SQCD from Seiberg-Witten curves

## I Introduction

The gaugino and meson vacuum condensates of SQCD are fundamental features of the theory, which are associated with the presence of a mass gap, confinement and instantons. Over the years, various methods have emerged for calculating the values of these condensates.

An indirect but simple way is to calculate the minima of one of two related non-perturbative quantum superpotentials: the Veneziano-Yankielovicz-Taylor (VYT) superpotential Taylor et al. (1983), which was conjectured on the basis of reproducing the correct $`U(1)`$ anomalies of the theory, or the Affleck-Dine-Seiberg (ADS) superpotential Affleck et al. (1984), which was deduced using holomorphy arguments, and the $`U(1)_A`$ selection rule. For the ADS superpotential, which only yields the meson condensate, one can obtain the gaugino condensate by using a suitable Konishi anomaly equation. Note that both of these superpotentials only exist when the number of flavours is less than the number of colours $`N_f<N_c`$. For the region $`N_fN_c`$, alternative techniques for calculating the condensates are needed.

There are also direct methods, as described in Konishi and Ricco (2003). The so-called strong coupling instanton (SCI) and weak coupling instanton (WCI) one-loop calculations Dorey et al. (2002), though developed for calculating the gaugino condensate in pure $`𝒩=1`$ SYM, can also be applied to SQCD by taking into account the differing one-loop beta function. The meson condensate can then be obtained from a Konishi anomaly equation. The results of the two methods disagree by a factor of 4/5, though it is widely believed that the SCI calculation is not reliable Dorey et al. (2002) and that we should believe the WCI calculation, which also agrees with the results from the VYT and ADS superpotentials. The main problem with the SCI approach seems to be that in the strong coupling confining phase of the theory, semi-classical instanton calculations are not valid. In another interesting technique Davies et al. (1999) which could be generalized to the case of fundamental matter, space-time is compactified on $`R^3\times S^1`$, in order that the semi-classical instanton calculation is valid.

In our approach, we use the intimate relationship between $`𝒩=1`$ and $`𝒩=2`$ SYM, extending Konishi and Ricco’s calculation Konishi and Ricco (2003) to the case of fundamental matter with $`SU(N_c)`$ gauge group and $`N_f<2N_c`$ flavours. In turn, their work was a generalization of an article by Finnell and Pouliot Finnell and Pouliot (1995) for $`SU(2)`$. It may come as no surprise that our results agree with the prediction of the VYT and ADS superpotentials, and the WCI calculation. Nevertheless, it is satisfying, to have confirmation from another independent source. Our method is also valid for a wider range of flavours than for the two superpotentials.

## II The Calculation

As previously mentioned, starting with $`𝒩=2`$ super Yang-Mills with $`SU(N_c)`$ gauge group and $`N_f<2N_c`$ flavours, we add a mass term $`\mu \text{Tr }\varphi ^2`$ for the $`𝒩=1`$ adjoint scalar field $`\varphi `$. In terms of the resultant $`𝒩=1`$ theory, the classical superpotential is given by

$$W_{tree}=\mu \text{Tr }\varphi ^2+\sqrt{2}\text{Tr }\underset{f=1}{\overset{N_f}{}}\stackrel{~}{Q}_f\varphi Q^f+\text{Tr }\underset{f=1}{\overset{N_f}{}}\stackrel{~}{Q}_fm_fQ^f,$$
(1)
where $`Q^f`$ and $`\stackrel{~}{Q}_f`$ are quark superfields in the fundamental representation of the $`SU(N_c)`$ gauge group, and $`m_f`$ are the quark masses for $`f=1,\mathrm{},N_f`$. We focus on the case $`N_f<2N_c`$ where the $`𝒩=2`$ theory is asymptotically free, hence Seiberg-Witten curves can be used.

In the limit $`\mu \mathrm{}`$, the field $`\varphi `$ decouples and can be integrated out to leave pure SQCD. Classically, when all $`N_f`$ quark masses are non-zero, the minima of the SQCD potential is given by $`Q^f=\stackrel{~}{Q}_f=0`$.

The Seiberg-Witten curve for pure $`𝒩=2`$ plus fundamental matter and $`SU(N_c)`$ gauge group is

$$y^2=\frac{1}{4}\underset{a=1}{\overset{N_c}{}}(x\varphi _a)^2\mathrm{\Lambda }^{2N_cN_f}\underset{f=1}{\overset{N_f}{}}(x+m_f),$$
(2)
where $`\varphi =\text{diag }(\varphi _1,\mathrm{}\varphi _{N_c})`$ and $`_{a=1}^{N_c}\varphi _a=0`$.

We can in principle calculate $`<Tr\varphi ^2>`$ for the $`𝒩=1`$ theory described by $`W_{tree}`$ in (1), using the Seiberg-Witten curve. The mass term $`\mu \text{Tr }\varphi ^2`$ lifts the vacuum degeneracy and causes the moduli, $`\varphi _a`$, to move to the point on the curve where the maximum number of massless magnetic monopoles condense. These $`(N_c1)`$ massless monopoles are responsible, via the Higgs mechanism, for giving a mass to the abelian gauge-fields belonging to the $`𝒩=2`$ multiplet associated with $`\varphi _a`$, as argued in Seiberg and Witten (1994a, b).

The maximum $`(N_c1)`$ massless magnetic monopoles correspond to $`(N_c1)`$ double zeros of the SW curve. So $`\varphi _a`$ are chosen such that the curve is of the form

$$y^2=F_2(x)H_{N_c1}^2(x),$$
(3)
where the subscript on $`F`$ and $`H`$ refers to the degree of the polynomial in $`x`$. The tuning of $`\varphi _a`$ such that the curve is of the above form is known as performing a complete factorization of the curve.

In practice, the curve with fundamental matter is extremely difficult to factorize. In fact, it is only recently that the general solution has been found Demasure and Janik (2003), by exploiting new powerful tools Dijkgraaf and Vafa (2002) for calculating F-terms of the $`𝒩=1`$ quantum superpotential. We will not need these techniques though, as in the limit $`\mu \mathrm{}`$, the curve (2) simplifies to a more manageable form, where the technique of Chebysev polynomials is applicable, as in Konishi and Ricco (2003).

First we write the $`𝒩=2`$ dynamical scale $`\mathrm{\Lambda }_{𝒩=2}`$ in terms of the scale on the curve $`\mathrm{\Lambda }`$

$$\mathrm{\Lambda }^2=2^{N_c/(N_c\frac{1}{2}N_f)}\mathrm{\Lambda }_{𝒩=2}^2.$$
(4)
We then match the dynamical scales of the two theories in the usual way

$$\mathrm{\Lambda }_{𝒩=1}^{3N_cN_f}=\mu ^{N_c}\mathrm{\Lambda }_{𝒩=2}^{2N_cN_f}.$$
(5)
Before taking the large $`\mu `$ limit, it is convenient to define new variables

$`\widehat{\varphi }_a=\sqrt{\mu }\varphi _a,`$ $`\widehat{x}=\sqrt{\mu }x,`$ $`\widehat{y}^2=\mu ^{N_c}y^2.`$ (6)
In the limit $`\mu \mathrm{}`$, while keeping $`\mathrm{\Lambda }_{𝒩=1}`$ constant, the Seiberg-Witten curve is given by

$$\widehat{y}^2=\frac{1}{4}\underset{a=1}{\overset{N_c}{}}(\widehat{x}\widehat{\varphi }_a)^2\frac{\mathrm{\Lambda }_{𝒩=1}^{3N_cN_f}}{2^{N_c}}detm,$$
(7)
where $`m=\text{diag }(m_1,\mathrm{}m_{N_f})`$, and $`\widehat{y}^2`$ can be factorized using Chebysev polynomials Douglas and Shenker (1995), just as if there’s no fundamental matter.

At the point where there is a maximum number of double zero’s, the curve becomes

$`\widehat{y}^2={\displaystyle \frac{\mathrm{\Lambda }_{𝒩=1}^{3N_cN_f}detm}{2^{N_c}}}(T_{N_c}^2(\xi )1),`$ (8)
where

$`\xi =({\displaystyle \frac{1}{2^{N_c}\mathrm{\Lambda }_{𝒩=1}^{3N_cN_f}detm}})^{\frac{1}{2N_c}}\widehat{x}e^{2\pi ik/2N_c}`$ $`k=1,\mathrm{},N_c,`$ (9)
where $`T_{N_c}(\xi )`$ is a Chebysev polynomial of the first kind, as described in appendix A, and we have used equation (20) to normalize $`\xi `$. So the rescaled moduli are

$$\widehat{\varphi }_a=\omega _a(2^{N_c}\mathrm{\Lambda }_{𝒩=1}^{3N_cN_f}detm)^{\frac{1}{2N_c}}e^{2\pi ik/2N_c},$$
(10)
where $`\omega _a`$ is defined in equation (17). Note that $`_a\widehat{\varphi }_a=0`$, since $`_a\omega _a=0`$. Furthermore

$$<\text{Tr }\widehat{\varphi }^2>=\underset{\mu \mathrm{}}{lim}\mu <\text{Tr }\varphi ^2>=N_c(\mathrm{\Lambda }_{𝒩=1}^{3N_cN_f}detm)^{\frac{1}{N_c}}e^{2\pi ik/𝒩_c},$$
(11)
where we have used that $`_{a=1}^{N_c}\omega _a\omega _a=N_c/2`$.

Using the Konishi anomaly equations Cachazo et al. (2003)

$`{\displaystyle \frac{1}{16\pi ^2}}<\text{Tr }\lambda \lambda >+<\text{Tr }q^f\varphi \stackrel{~}{q}_f>`$ $`=`$ $`{\displaystyle \frac{\mu }{N_c}}<\text{Tr }\varphi ^2>,`$ (12)
$`<\text{Tr }q^f\stackrel{~}{q}_g>+\sqrt{2}(m^1)_{}^{f}{}_{g}{}^{}<\text{Tr }\stackrel{~}{q}_g\varphi q^f>`$ $`=`$ $`{\displaystyle \frac{(m^1)_{}^{f}{}_{g}{}^{}}{16\pi ^2}}<\text{Tr }\lambda \lambda >,`$ (13)
where $`\lambda `$ is the gaugino field, and $`q^f`$ and $`\stackrel{~}{q}_g`$ the lowest components of the quark superfields $`Q^f`$ and $`\stackrel{~}{Q}_g`$, we calculate the gaugino and meson condensates. Note that in the limit $`\mu \mathrm{}`$, using equation (5), $`\mathrm{\Lambda }_{𝒩=2}`$ becomes zero and hence so does $`\varphi _a`$. Thus, terms of the form $`<\text{Tr }q\varphi \stackrel{~}{q}>`$ can be ignored.

Finally, the gaugino and meson vacuum condensates are given by

$`{\displaystyle \frac{1}{16\pi ^2}}<\text{Tr }\lambda \lambda >=(\mathrm{\Lambda }_{𝒩=1}^{3N_cN_f}detm)^{\frac{1}{N_c}}e^{2\pi ik/𝒩_c},`$ (14)
$`<\text{Tr }q^f\stackrel{~}{q}_g>=(m^1)_{}^{f}{}_{g}{}^{}(\mathrm{\Lambda }_{𝒩=1}^{3N_cN_f}detm)^{\frac{1}{N_c}}e^{2\pi ik/𝒩_c},`$ (15)
respectively, where $`k=1\mathrm{}N_c`$.

###### Acknowledgements.

I would like to thank Paolo Di Vecchia for suggesting the problem, and Yves Demasure, Paolo Merlatti and Gabriele Ferretti for helpful discussions. This work was partly funded by a Marie Curie training site fellowship.

## Appendix A The Chebysev polynomials

The Chebysev polynomials of the first and second kind are $`T_n(\xi )=\mathrm{cos}(n\mathrm{cos}^1\xi )`$ and $`U_n(\xi )=[\mathrm{sin}((n+1)\mathrm{cos}^1\xi )/\sqrt{1\xi ^2}]`$ respectively. They obey the relation

$$T_n^2(\xi )1=(\xi ^21)U_{n1}^2(\xi )$$
(16)
The zero’s of $`T_n`$ are given by

$`\omega _a=\mathrm{cos}({\displaystyle \frac{\pi }{n}}(a{\displaystyle \frac{1}{2}}))`$ $`a=1,\mathrm{},n,`$ (17)
and of $`U_n`$ by

$`\zeta _a=\mathrm{cos}({\displaystyle \frac{\pi }{n+1}}a)`$ $`a=1,\mathrm{},n,`$ (18)
Thus,

$`{\displaystyle \underset{a=1}{\overset{n}{}}}\omega _a=0,`$ $`{\displaystyle \underset{a=1}{\overset{n}{}}}\omega _a\omega _a=n/2.`$ (19)
It is useful to know the explicit factorized form of $`T_n`$ and $`U_n`$

$$T_n(\xi )=\{\begin{array}{cc}\frac{1}{2}_{a=1}^{n/2}(4\xi ^24\omega _a^2)\hfill & \text{(n even)}\hfill \\ \frac{1}{2}(2\xi )_{a=1}^{(n1)/2}(4\xi ^24\omega _a^2)\hfill & \text{(n odd)},\hfill \end{array}$$
(20)
$$U_n(\xi )=\{\begin{array}{cc}_{a=1}^{n/2}(4\xi ^24\zeta _a^2)\hfill & \text{(n even)}\hfill \\ (2\xi )_{a=1}^{(n1)/2}(4\xi ^24\zeta _a^2)\hfill & \text{(n odd)}.\hfill \end{array}$$
(21)


text:
# Blue horizontal branch stars in the Sloan Digital Sky Survey: I. Sample selection and structure in the Galactic halo

## 1 Introduction

Blue horizontal branch (BHB) stars have long been important objects of study as tracers of the Galactic potential (Pier, 1983; Sommer-Larsen, Christensen, & Carter, 1989; Clewley et al., 2002). Their constant absolute magnitude allows accurate space positions to be determined. Historically, many researchers have used samples of BHB stars, or other standard-candle stars such as RR Lyraes, to study the local halo, i.e. the vicinity of the solar neighborhood where typical distances are of order the solar Galactocentric radius $`R_0=8\mathrm{kpc}`$. These samples can be used to constrain the kinematics of the local halo (Pier, 1984; Layden et al., 1996). The distant halo (i.e. typical distances $`R_0`$), however, is less well constrained (but see Sommer-Larsen et al. (1997)). Isolating a sample of distant BHB stars would enable study of the structure and kinematics of the Galactic halo on larger scales (Sirko et al., 2003, hereafter Paper II), and estimates of the mass of the Galaxy.

The Sloan Digital Sky Survey (SDSS), currently halfway through a photometric and spectroscopic survey of a quarter of the sky, offers an excellent opportunity to isolate and study BHBs. The SDSS spectroscopic survey targets objects between about $`g`$ = 15.5 to 19 (the bright limit is set to avoid contamination in adjacent spectroscopic fibers at the slit end of the fiber bundle), although a significant number of stars outside this range is also observed. Thus BHB stars at distances from $`7`$ to $`60\mathrm{kpc}`$ from the sun appear in the spectroscopic survey.

This paper presents $`1170`$ BHB stars selected from SDSS that probe the outer halo. The plan of the paper is as follows. Pertinent characteristics of the Sloan Digital Sky Survey are summarized in §2. §3 describes the preliminary color cut. §4 discusses spectroscopic BHB classification criteria. For the faint stars, the spectroscopic classification breaks down so we resort to a more stringent color cut described in §5. Radial velocity errors are evaluated in §6. §7 discusses the calculation of bolometric corrections for accurate distances, and evaluates magnitude errors. Some sanity checks are presented in §8, including an analysis of the limited proper motion data available, and a map of the Galactic halo which shows the conspicuous Sagittarius tidal stream (figure 13). The results are summarized in §9. The sample is available in Table 3 in the electronic version of this paper.

## 2 The Sloan Digital Sky Survey

The Sloan Digital Sky Survey (SDSS: York et al., 2000) is a project to image about 1/4 of the sky at high Galactic latitude to a depth of about 21 to 22 magnitudes in five broad bands covering the entire optical range and to obtain spectra of about a million galaxies and 100,000 quasars selected from the imaging data. While galaxies and quasars are the primary targets of the spectroscopic survey, the SDSS also obtains a vast amount of stellar data. Stellar spectra are acquired by several routes: observations of calibration stars, observations of various classes of star to be used as “filler” observations when a given region of sky has too low a density of primary targets, and via objects targeted as candidate galaxies and quasars that turn out to be stars.

The SDSS uses a dedicated 2.5 meter telescope and a large format CCD camera (Gunn, Carr, Rockosi, & Sekiguchi et al., 1998) at the Apache Point Observatory in New Mexico to obtain images almost simultaneously in five broad bands ($`u,g,r,i,z`$) centered at 3551, 4686, 6166, 7480, and 8932 Å respectively (Fukugita et al., 1996; Gunn, Carr, Rockosi, & Sekiguchi et al., 1998). The imaging data are automatically processed through a series of software pipelines which find and measure objects and provide photometric and astrometric calibrations to produce a catalogue of objects with calibrated magnitudes, positions and structure information. The photometric pipeline Photo (Lupton et al., 2001, and in preparation) detects the objects, matches the data from the five filters, and measures instrumental fluxes, positions and shape parameters. The last allow the classification of objects as “point source” (compatible with the point spread function \[psf\]), or “extended.” The instrumental fluxes are calibrated via a network of primary and secondary stellar flux standards to $`\mathrm{AB}_\nu `$ magnitudes (Fukugita et al., 1996; Hogg et al., 2001; Smith et al., 2002). We use the “psf magnitude,” which represents the magnitude of the best-fit psf to the object, as determined by Photo version 5.4.25 (2003 June). The psf magnitudes are currently accurate to about 2% in $`g,r`$, and $`i`$ and 3-5% in $`u`$ and $`z`$ for bright ($`<`$ 20 mag) point sources. The SDSS is 50% complete for point sources to $`(u,g,r,i,z)=(22.5,23.2,22.6,21.9,20.8)`$, and the full-width at half-maximum (FWHM) of the psf is about $`1.5^{\prime \prime }`$ (Abazajian et al., 2003). The data are saturated at about 14 mag in $`g`$, $`r`$ and $`i`$ and about 12 mag in $`u`$ and $`z`$. Astrometric calibration is carried out using a second set of less sensitive CCDs in the camera, which allows the transfer of astrometric catalogue positions to the fainter objects detected by the main camera. Absolute positions are accurate to better than $`0.1^{\prime \prime }`$ in each coordinate (Pier et al., 2003).

The SDSS spectroscopic survey is carried out by targeting objects selected from the imaging catalogues according to pre-determined and carefully tested criteria (Stoughton et al., 2002). The targeted stars, quasars and galaxies have limiting magnitudes of around $`1820`$ mag (Richards et al., 2002; Eisenstein et al., 2001; Strauss et al., 2002). BHB stars appear in the SDSS spectroscopic sample via two main routes. At redshifts around 2.7, the quasar color locus crosses the BHB color locus (Richards et al., 2002) and many BHB stars are targeted as candidate quasars. BHB stars are also photometrically selected to serve as filler observations in regions of the sky with a low density of primary targets. The sky locations of the target objects are mapped to a series of plug plates which feed a pair of double spectrographs (A. Uomoto, S. Smee et al. unpublished) via optical fibers. The spectrographs jointly observe 640 objects (including sky positions and standard stars) simultaneously, and provide continuous wavelength coverage from about 3800 Å to 9500 Å at a resolution of about $`18002100`$. The nominal total exposure time for each spectroscopic plate is 45 minutes (15 minutes $`\times `$ 3 exposures), plus observations of flat field and arc standard lamps for every plate. The entrance aperture of the fibers is $`3^{\prime \prime }`$. The spectra are matched to the photometric objects using a video mapping system for the plug plates (D. Schlegel and D. Finkbeiner, unpublished). The spectra are recorded on CCD detectors. The data are optimally extracted, the sky spectrum is subtracted, and wavelength, flat-field and spectrophotometric calibrations are applied automatically by the Spectro-2D software pipeline, written by D. Schlegel and S. Burles (unpublished).

Each spectroscopic plate is a circle of diameter $`3^{}`$ on the sky. The present work is based on the data of $`937`$ unique spectroscopic plates taken to date, through 2003 Aug 6. In addition, $`77`$ plates were observed redundantly, but in defining the BHB sample we exclude redundant observations by adopting the plate with the highest S/N. However, we use the information from these redundant observations in §6 below to evaluate velocity errors.

The radial velocity is found by fitting the spectrum to template spectra, using an automated fitting program (Spectro\_Brightstars) written by D. Schlegel (unpublished). The basic procedure is as follows. Prugniel & Soubiran (2001) have published a catalogue of about a thousand spectra of stars with a wide range of properties (spectral types, metallicities etc.) observed at high spectral resolution with the ELODIE spectrograph. One version of the ELODIE spectra is binned to R $``$ 20,000 and spectrophotometrically calibrated to about 5%. Brightstars was used to select template spectra from the SDSS data which match the ELODIE spectra. These template spectra are typically composites from many SDSS spectra, and are set up to provide relatively coarse sampling of the spectral type grid. The template spectra are aligned to zero heliocentric velocity. Each SDSS spectrum is then compared with the template spectra in wavelength space, and the template and wavelength shift that give the best match is automatically measured. The stars in the present sample were almost all best-fit with spectral types between B6 and A0. The radial velocities are corrected to the heliocentric standard of rest, and throughout this paper are denoted $`v_r`$.

## 3 Color cut

BHB stars have very characteristic colors in the high-Galactic latitude stellar sample observed by the SDSS imaging survey, and photometric selection of candidate objects is straightforward. BHB stars are bluer in $`gr`$ than most halo stars because the halo main sequence turnoff lies at spectral type F/G. The large Balmer jump in BHB low-gravity stars gives them redder $`ug`$ colors than other blue objects such as low redshift quasars and white dwarfs, and indeed the SDSS $`u`$ and $`g`$ filters are specifically designed to provide optimal photometric separation of low-redshift quasars (Gunn, Carr, Rockosi, & Sekiguchi et al., 1998; York et al., 2000; Richards et al., 2002). The SDSS $`gr`$ vs. $`ug`$ color-color diagram in figure 1 shows all objects in the spectroscopic sample of $`937`$ plates used herein that were subsequently identified as stars by Brightstars, which is accurate enough that quasars and galaxies are virtually nonexistent in the sample. (Throughout this paper all magnitudes are corrected for Galactic extinction using the maps of Schlegel, Finkbeiner, & Davis (1998).) The region occupied by BHB stars is indicated with a “color cut box,” with boundaries $`0.8<ug<1.35,0.4<gr<0.0`$. The boundaries of this selection box are similar to those used by Yanny & Newberg et al. (2000), and figure 8 (see §5) shows that our final sample of bright BHB stars is well contained within the boundaries. Stars that lie in this region of color-color space pass what we will refer to as the “color cut.”

Although BHB stars are effectively selected through the color cut described above, there remains a significant degree of contamination from two other types of star. The locus of cooler F stars is actually quite outside of the color cut box, to the red in $`gr`$ (Yanny & Newberg et al., 2000), but there are so many F stars that they leak into the BHB candidate sample, via intrinsic variations and photometric errors, with numbers comparable to the number of BHB stars. To identify these cooler stars merely requires a measure of temperature independent of $`ug`$ or $`gr`$. Somewhat more challenging to identify are the A-type stars with higher surface gravity (Preston & Sneden, 2000; Yanny & Newberg et al., 2000; Clewley et al., 2002). These are probably field blue stragglers, which are some 2 magnitudes less luminous (Preston & Sneden, 2000). While the Balmer jump as measured by $`ug`$ can separate out many of these stars, some remain in the sample. They can be removed using their broader Balmer lines as described in the next section.

## 4 Balmer line analysis cuts

The upper panel of figure 2 presents a typical BHB spectrum, showing the deep Balmer lines of an A-type star. In this section we describe two methods for analyzing these lines: the widely used $`D_{0.2}`$ method (Pier, 1983; Sommer-Larsen & Christensen, 1986; Arnold & Gilmore, 1992; Flynn, Sommer-Larsen, & Christensen, 1994; Kinman, Suntzeff, & Kraft, 1994; Wilhelm, Beers, & Gray, 1999) and the scalewidth-shape method (Clewley et al., 2002). We focus on the H$`\gamma `$ and H$`\delta `$ lines because the the H$`\beta `$ and H$`\alpha `$ lines have lower continuum levels, and the lines above H$`\delta `$ are much more closely spaced, making the determination of the continuum difficult (Yanny & Newberg et al., 2000).

To analyze the lines, the continuum must first be fitted and divided out. For the $`D_{0.2}`$ method, following Yanny & Newberg et al. (2000), we take the continuum to range from 4000 to $`4500\text{Å}`$ with $`60\text{Å}`$ masks covering the positions of the H$`\gamma `$ and H$`\delta `$ lines. We then divide the spectrum by the best-fit sixth-order Legendre polynomial to this extracted continuum. Finally, the spectrum is boxcar-smoothed with a boxcar length of 5 pixels. An example normalized and extracted SED is shown in the lower panel of figure 2. Similarly, for the scalewidth-shape method, we take the continuum extraction rule of Clewley et al. (2002) (which ranges from 3863 to $`4494\text{Å}`$) and fit to a fourth-order Legendre polynomial. After the continuum is divided out with this prescription, the SED looks very similar to the example SED in the lower panel of figure 2.

The following terminology is used in this section and throughout. A “method” refers to a general algorithm, which may or may not be desirable for isolated use, which employs one or more “criteria” for discriminating BHB stars. The $`D_{0.2}`$ method and the scalewidth-shape method are presented in the following. A “cut” is a method or combination of methods that we *actually* use to eliminate contaminants from the sample. The combination cut discussed below incorporates the $`D_{0.2}`$ and scalewidth-shape methods.

### 4.1 The $`D_{0.2}`$ method

Since BHB stars have low surface gravity, their Balmer lines are narrower than the smaller main sequence blue stragglers. The $`D_{0.2}`$ method discriminates BHB stars from BS stars by determining the value of $`D_{0.2}`$, the width in angstroms of the Balmer line at 80% of the continuum Yanny & Newberg et al. (2000). Figure 3 demonstrates how this parameter can discriminate a BHB from a BS star. Although simple, this method has been found to be fairly robust and can provide a better measurement of the linewidth than other more complicated procedures (e.g. Flynn, Sommer-Larsen, & Christensen, 1994).

The $`D_{0.2}`$ method nominally only applies to objects with A-type spectra. However, as mentioned in $`\mathrm{\S }`$3, a large contaminant of our sample is F stars. Since F-type stars have weaker Balmer lines than A-type stars, one way of measuring temperature is by measuring the depth of the H$`\gamma `$ or H$`\delta `$ line. Here we define $`f_\mathrm{m}`$ as the flux relative to the continuum at the minimum of the line. Thus, higher values of $`f_\mathrm{m}`$ for the Balmer lines signify cooler stars. A plot of $`D_{0.2}`$ vs. $`f_\mathrm{m}`$ for stars passing the color cut is shown in figure 4. In this figure, only stars with $`g`$ magnitudes $`<18`$ are shown, so that the trends in the plot are not obscured by the noise in the data of fainter stars (see figure 5, to be discussed in §4.3).

The concentration of stars near $`(f_\mathrm{m},D_{0.2})=(.30,26\text{Å})`$ represents the BHB stars, and those stars with larger values of $`D_{0.2}`$ are BS stars. Figure 4 also suggests that the locii of BS and BHB stars propagate to larger $`f_\mathrm{m}`$ and smaller $`D_{0.2}`$ until the two locii converge at about $`(f_\mathrm{m},D_{0.2})=(.4,23\text{Å})`$. In other words, for a given $`.3f_\mathrm{m}.4`$, the distribution of stars in $`D_{0.2}`$ is bimodal, but the stars are most cleanly separated for the smaller values of $`f_\mathrm{m}`$ and not separated at all when $`f_\mathrm{m}.4`$. Even though figure 4 doesn’t show the trend clearly, we return to this point in §4.3 and show that the stars with smaller values of $`D_{0.2}`$ but with $`.34f_\mathrm{m}.4`$ are the other horizontal branch stars; RR Lyraes, for example, fall into this region of the horizontal branch. Because horizontal branch stars are evidently more difficult to discriminate from main sequence stars at lower temperatures, we adopt the criteria $`f_\mathrm{m}<.34,D_{0.2}<28`$ for likely BHB stars for the H$`\delta `$ $`D_{0.2}`$ method.

### 4.2 The scalewidth-shape method

Clewley et al. (2002) have proposed the scalewidth-shape method for discriminating BHB stars based on the Sérsic profile for Balmer lines (Sérsic, 1968):

$$y=1.0a\mathrm{exp}\left[\left(\frac{|xx_0|}{b}\right)^c\right].$$
(1)
Here, $`y`$ is the normalized flux and $`x`$ is the wavelength in Å. We fit the normalized extracted spectrum (where the bounds of extraction are given by Clewley et al. (2002)) to the Sérsic profile with three free parameters: $`a,b`$ and $`c`$. The parameter $`x_0`$ is assumed to be the nominal location of the Balmer line corrected for redshift, which has already been determined. After fitting both Balmer lines of every star to the three-parameter ($`a,b,c`$) Sérsic profile, there is a reasonably small dispersion in the determination of $`a`$: for the $`H\gamma `$ line the average is $`.690`$ and the root-mean-square (rms) is $`.019`$; for $`H\delta `$ the average is $`.740`$ and the rms is $`.018`$. We fix $`a`$ at its average value among stars passing the BHB criteria of the H$`\delta `$ $`D_{0.2}`$ method, which is not objectionable because the $`D_{0.2}`$ method is independent of the scalewidth-shape method. The Balmer lines are then refitted to the two-parameter ($`b,c`$) Sérsic profile. This refitting is important because errors in $`a`$ are correlated with the errors in $`b`$ and $`c`$ (Clewley et al., 2002). For this work we adopt $`a_\gamma =.690`$ for the H$`\gamma `$ line and $`a_\delta =.740`$ for the H$`\delta `$ line. The optimum value of $`a`$ most likely depends on the resolution of spectroscopy and similar factors.

Fitting to the Sérsic profile thus provides the two parameters $`b`$ and $`c`$ for both the H$`\gamma `$ and H$`\delta `$ lines, for every star. Clewley et al. (2002) show that BHB and BS stars separate rather cleanly in a plot of $`b`$ against $`c`$. Furthermore, they show that $`c`$ is a measure of temperature: cooler stars have smaller values of $`c`$.

### 4.3 Determination of criteria for the $`D_{0.2}`$ and scalewidth-shape methods

We concentrate on the set of brighter stars with $`g<18`$ because the noise in the spectra of fainter stars throws large errors in the parameter determinations for the $`D_{0.2}`$ method and the scalewidth-shape method. This point is illustrated in figure 5. We will return to the selection of BHBs from the fainter star sample in §5.

Figure 6 shows the parameters of the sample of BHB candidates with $`g<18`$ as determined by the two linewidth analysis methods described in the preceding sections. Each star is represented by the same color in all four panels. In this figure the H$`\delta `$ $`D_{0.2}`$ method is represented in the lower left panel, and stars which pass the BHB criteria for this method are colored blue, stars which are likely blue stragglers are green, and all others (cooler stars) are red. It is readily apparent that the four *independent* parameter spaces map BHB stars, BS stars, and cool stars to different regions.

The exact BHB criteria for each method were determined through a combination of adjusting the bounds of the criteria by eye and comparing the results of one method to the others. For instance, the scalewidth-shape method separates BHB and BS stars cleanly enough that there is a natural gap between the two distributions in a plot of $`b`$ vs. $`c`$, as shown in figure 6. A parabola is qualitatively drawn to divide the two distributions. Similarly, with regard to the $`D_{0.2}`$ method, the distribution in $`D_{0.2}`$ of the stars with low values of $`f_\mathrm{m}`$ naturally separates into the two populations of BHB and BS stars. The issue of where to draw the temperature criterion ($`f_\mathrm{m}`$ for the $`D_{0.2}`$ method and $`c`$ for the scalewidth-shape method) is more sensitive because the “natural” gap is harder to see in the $`D_{0.2}`$ method and seems to be nonexistent in the scalewidth-shape method. In fact, the upper right panel of figure 6 hints that the division between stars with low surface gravity and high surface gravity continues even for much cooler stars than BHB stars. The cooler population ($`0.8c_\gamma 1.0`$) with lower surface gravity (lower $`b_\gamma `$) is thus the horizontal branch stars which are not blue, e.g. RR Lyraes. This $`b`$ vs. $`c`$ plot corroborates the evidence in §4.1 that non-blue horizontal branch stars form a locus in the $`(f_\mathrm{m},D_{0.2})`$ plane which extends out to cooler temperatures corresponding to $`f_\mathrm{m}.4`$. Figure 7 shows that the stars of the population that has smaller values of $`b_\gamma `$ also have smaller values of $`D_{0.2}`$. With this information, we deduce that the BHB stars are located in the clump in the H$`\delta `$ $`f_\mathrm{m}D_{0.2}`$ plane centered on $`(f_\mathrm{m},D_{0.2})_{\mathrm{H}\delta }(.30,26\text{Å})`$ and we set the temperature criterion at $`f_\mathrm{m}=.34`$. A subsequent inspection of figure 6 reveals that most BHB stars have $`c_\gamma >1.0`$, so we adopt this as the temperature criterion for the H$`\gamma `$ scalewidth-shape method. Table 1 lists the adopted BHB selection criteria for each method.

### 4.4 The combination cut

We define the combination cut as those stars which pass the criteria for (the H$`\gamma `$ $`D_{0.2}`$ and H$`\delta `$ scalewidth-shape methods) OR (the H$`\delta `$ $`D_{0.2}`$ and H$`\gamma `$ scalewidth-shape methods). We assume that the H$`\gamma `$ $`D_{0.2}`$ and H$`\delta `$ scalewidth-shape methods are completely independent of each other, as are the H$`\delta `$ $`D_{0.2}`$ and H$`\gamma `$ scalewidth-shape methods. This assumption is reasonable because no information is shared in the analyses of either line, and the $`D_{0.2}`$ and scalewidth-shape methods measure fundamentally different properties of the spectral energy distribution. Under this assumption, the contamination fraction for the combination cut would be $`\kappa _{\mathrm{combo}}2\kappa ^2`$, where $`\kappa `$ is the contamination fraction of a single method. In Appendix A we show that the contamination fraction of any one of the four methods is probably $`\kappa 10\%`$, so $`\kappa _{\mathrm{combo}}2\%`$ under these assumptions. However, to be more conservative, the fact that $`\kappa 10\%`$ for all four individual methods suggests that $`\kappa _{\mathrm{combo}}`$ is probably $`10\%`$ or less; in Paper II we use $`\kappa _{\mathrm{combo}}=10\%`$ to study errors introduced by contamination into the kinematic analysis .

## 5 The fainter stars

As shown in figure 5, the linewidth analysis methods cannot be expected to perform accurately for fainter stars than $`g18`$; in fact the number of $`g>18`$ stars that are identified differently in any two methods ($`100400`$) is of the same order as the number of $`g>18`$ stars identified as BHB stars in one method ($`200500`$). Because these fainter stars have noisier spectra, not only are the algorithmic BHB criteria unreliable, but eye-inspection would be difficult at best. However, an inspection of the results of the combination cut applied to the bright ($`g<18`$) stars, as shown in figure 8, reveals that BHB stars, BS stars, and cooler stars occupy different regions of color-color space, although there is significant overlap. In this figure BHB stars are defined as those which pass the combination cut; (green) BS stars are defined as those which do not pass the combination cut but lie in the BS region for the H$`\delta `$ $`D_{0.2}`$ method; and red stars are everything else. The thick curved line with periodic circles will be explained in §7.

The bounds of this particular plot are the same as the color cut discussed in §3, so it is evident that the BHB stars are more likely to be found in the region enclosed by the “piano” shape indicated than outside of it. Since SDSS photometry can be expected to be accurate to a fainter magnitude ($`g20`$, York et al., 2000) than can the combination method (see figure 5), we select BHB stars with $`g>18`$ solely on the basis of color. The new color cut is indicated by the piano shape in figure 8 (and is given by the region $`0.85<ug<1.3,0.31<gr<.13`$ excluding the elliptical subregion $`((ug0.85)/0.31)^2+((gr+0.13)/0.11)^2<1`$) and will be referred to here as the “stringent” color cut.

What are the contamination and incompleteness for the stringent color cut with no spectroscopic analysis? Assuming that the combination cut is 100% accurate for $`g<18`$, there are $`733`$ BHB stars, $`172`$ of which do not pass the stringent color cut. There are $`199`$ stars that pass the stringent color cut that are not classified as BHB stars by the combination method. Most of these ($`133`$) are classified as blue stragglers by the H$`\delta `$ $`D_{0.2}`$ method. Thus the contamination for the stringent color cut alone is $`199/(733172+199)25\%`$.

## 6 Duplicate stars and radial velocity errors

Since the spectroscopic plates overlap on the sky, sometimes objects are targeted for spectroscopy on two different plates. Of the $`4515`$ stars passing the color cut (from *unique* plates), $`81`$ are duplicate observations of the same star on different plates. We remove these $`81`$ redundant observations from our BHB sample. In addition, from the data of the $`77`$ non-unique plates (see §2), there are $`410`$ stars passing the color cut which duplicate observations of stars from the unique plates. Each spectroscopic observation of a star provides an independent measurement of radial velocity $`v_r`$ and, for the bright ($`g<18`$) stars, an independent classification as a BHB star or non-BHB star. Figure 9 shows the $`491`$ stars which have multiple spectra and therefore independent determinations of $`v_r`$. The ordinate axis shows the difference in the radial velocity determinations $`\mathrm{\Delta }v_r`$ between each duplicate star and its corresponding star in our purified sample of unique stars. Stars which are classified as BHBs for both observations are colored blue, stars which are classified as non-BHBs for both observations are colored red, and stars which have conflicting classifications between the two observations are colored green. Brighter than $`g=18`$, there are $`91`$ stars which are doubly classified as BHBs, $`131`$ stars which are doubly classified as non-BHBs, and $`13`$ stars which have conflicting classifications. This is consistent with the $`10\%`$ contamination determined in Appendix A, assuming that there is about an equal number of contaminants in the BHB sample as the number of BHB stars which are not classified as such (i.e., assuming $`\kappa \eta `$). Figure 9 also shows the dependence of velocity errors on magnitude; as expected, fainter stars tend to have greater velocity errors, presumably due to noiser spectra. The rank correlation coefficient of $`(g,|\mathrm{\Delta }v_r|)`$ is $`0.263`$. The velocity differences determined this way turn out to be consistent with the errors derived from the spectral fits for only about half of the stars, independent of the magnitude or whether or not a star is a BHB or BS star. This comparison implies that the errors have been underestimated and this is currently under investigation. The rms of $`\mathrm{\Delta }v_r`$ in figure 9 is $`36.5\mathrm{km}\mathrm{s}^1`$. Clearly the errors in radial velocity measurements are not gaussian, but a meaningful approximate result may still be obtained if they are assumed to be so. Since the distribution of $`\mathrm{\Delta }v_r`$ is based on the errors of two measurements, the error of one measurement is approximately $`\mathrm{\Delta }v_r36.5/\sqrt{2}26\mathrm{km}\mathrm{s}^1`$.

## 7 The absolute magnitude of BHB stars

In order to refine our photometric distance estimates, we have computed theoretical colors and bolometric corrections for BHBs. We integrate the latest Kurucz model atmospheres (downloaded from www.kurucz.harvard.edu) against the SDSS filter response functions at a standard airmass of 1.2 atmospheres, in the manner described by Lenz et al. (1998). Unfortunately, the latter authors report only colors, so it was necessary to repeat their calculations in order to obtain the bolometric corrections. Apart from wavelength, the Kurucz spectra are functions of three parameters: effective temperature $`T_{\mathrm{eff}}`$, surface gravity $`g=GM/r^2`$, and metallicity relative to solar $`[M/H]`$.

To derive absolute magnitudes and bolometric corrections, it is necessary to assume a relationship between the luminosity $`L`$ and $`T_{\mathrm{eff}}`$. The models of Dorman, Rood, & O’Connell (1993) for the theoretical zero-age horizontal branch (ZAHB) at $`[M/H]1.5`$ can be summarized roughly by

$$\mathrm{log}(L/L_{})1.580.73\mathrm{log}(T_{\mathrm{eff}}/10^4\mathrm{K})$$
(2)
over the range $`8000\mathrm{K}T_{\mathrm{eff}}12000\mathrm{K}`$. From a detailed study of the horizontal branch of M5, however, Baev, Markov, & Spassova (2001) find that luminosity varies more weakly with $`T_{\mathrm{eff}}`$; in fact their results appear to be consistent with constant luminosity and mass,

$$M0.6M_{},\mathrm{log}\frac{L}{L_{}}1.65$$
(3)
with a scatter of no more than $`\pm 0.05M_{}`$ and $`\pm 0.05\text{dex}`$, respectively, over the same temperature range. We adopt the mean values in equation (3), so that $`\mathrm{log}g3.52+4\mathrm{log}(T_{\mathrm{eff}}/10^4\mathrm{K})`$.

To give an idea of the theoretical uncertainties, we note that if in fact the relation (2) holds for our HB stars, then by adopting (3) we have overestimated their distance moduli by an average of $`0.18\mathrm{mag}`$ over the indicated range of $`T_{\mathrm{eff}}`$, and by neglecting the variation of $`L`$ with $`T_{\mathrm{eff}}`$, we have increased the scatter in distance modulus by $`0.09\mathrm{mag}`$, which should be taken in quadrature with the photometric measurement error.

The results of this exercise are shown in Table 2. Evidently, for the colors of interest to us, $`T_{\mathrm{eff}}(8000,12000)\mathrm{K}`$, and the influence of metallicity is slight if $`[M/H]1`$.

The sequence of theoretical horizontal branch stars, parametrized by temperature, is projected onto the $`(ug,gr)`$ plane in figure 8. The five black circles on the thick line denote, starting at the top and progressing down and to the left, stars with absolute magnitudes $`g=(0.6,0.55,0.6,0.7,0.8)`$. For every observed BHB star in our sample (colored blue in the figure) we define the most probable absolute magnitude associated with its $`(ug,gr)`$ colors by simply finding the absolute magnitude of the point on the theoretical track that is closest to the observed star in this color-color space.

The stringent color cut (piano shape in figure 8, see §5) was defined before the theoretical track was calculated, so the fact that the theoretical track satisfies the stringent color cut is a meaningful consistency check.

As a further consistency check, we have isolated several Galactic globular clusters observed by the SDSS photometric survey to date and plotted their color-magnitude diagrams (CMDs). Pal 3, NGC 5904, and NGC 6205 do not show a discernable horizontal branch, so they are not considered here. The $`g`$ vs. $`gr`$ CMDs of Pal 5 (Yanny & Newberg et al., 2000; Ivezić et al., 2000) and NGC 2419 are presented in figure 10, where the distance modulus, derived from cluster distances of $`d=84.2`$ and $`23.2\mathrm{kpc}`$ for NGC 2419 and Pal 5, respectively (Harris, 1996), has been subtracted from the dereddened apparent magnitude $`g`$. Stars within $`8^{}`$ of the center of NGC 2419 and within $`14^{}`$ of the center of Pal 5 were selected. The SDSS automated photometry avoids analyzing very crowded fields, so the stars presented here, though within the stated radii, tend to lie in the outskirts of the clusters. Vertical lines at $`gr=(0.4,0.0)`$ indicate our color cut. The theoretical horizontal branch sequence is also plotted in both panels. The horizontal branch stars of NGC 2419 appear to be systematically fainter than the theoretical prediction, while those of Pal 5 seem to be systematically brighter; however, these systematic differences could easily be due to inaccurate distances $`d`$. By selecting from these two globular clusters the stars passing the color cut ($`0.8<ug<1.35,0.4<gr<0.0`$) and lying in the region of the CMD for the horizontal branch ($`0.0<g5\mathrm{log}(d/10\mathrm{pc})<1.5`$), we compare the directly observed absolute magnitude $`g5\mathrm{log}(d/10\mathrm{pc})`$ to the indirectly obtained theoretical absolute magnitude based on the $`(ug,gr)`$ colors. For NGC 2419, the rms of the differences of the two among the 94 stars is .134 magnitudes, and the root-sum-square (rss) is .167 magnitudes<sup>1</sup><sup>1</sup>1 $`\mathrm{rms}=\sqrt{\frac{1}{N}\mathrm{\Sigma }\left(\mathrm{\Delta }g_i\mathrm{\Delta }g\right)^2}`$, $`\mathrm{\Delta }g=\frac{\mathrm{\Sigma }\mathrm{\Delta }g_i}{N}`$, $`\mathrm{rss}=\sqrt{\frac{1}{N}\mathrm{\Sigma }\mathrm{\Delta }g_i^2}`$. For Pal 5, the rms is .126 magnitudes and the rss is .202 magnitudes. The difference between the rms and rss is due to the nonzero mean, or systematic error as discussed above. Therefore, the magnitude error of a star, with absolute magnitude determined by comparison to the theoretical horizontal branch sequence described above, is of order .2 magnitudes. The error may be as small as $`.13`$ magnitudes if the only cause of the systematic difference between directly observed absolute magnitudes and our indirectly obtained absolute magnitudes is incorrect distances for Pal 5 and NGC 2419.

## 8 Sanity checks

As figure 8 shows, BHB stars determined by the combination cut are in general redder in $`ug`$ and bluer in $`gr`$ than other objects such as blue stragglers. This result is consistent with the work of Yanny & Newberg et al. (2000), so it is reassuring. In this section we present several more sanity checks which verify that our sample of BHB stars is probably accurate.

### 8.1 Metallicity from the CaII K line

BHB stars are halo stars, and as such they should have lower metallicity than the disk star contaminants. The most promising spectral feature to use for a metallicity determination is the CaII K line at $`3933.7\text{Å}`$ (the H line, at $`3968.5\text{Å}`$, is not used because it is too close to H$`ϵ`$). Unfortunately, it is impractical to use the CaII K line as a *reliable* indicator of metallicity in individual SDSS spectra because much higher spectral resolution would be necessary. However, the following exercise is executed anyway as a sanity check. By employing a continuum-subtraction algorithm very similar to that described in §4.1 for the H$`\gamma `$ and H$`\delta `$ Balmer lines, we find the value of $`f_\mathrm{m}`$ for the CaII K line; smaller values indicate larger metallicity. Figure 11 presents the distributions of this parameter for our sample of BHB stars and non-BHB stars passing the color cut. The figure shows different distributions for bright ($`g<18`$) and faint ($`g>18`$) stars because of the difference in S/N. Reassuringly, our BHB stars tend to have lower metallicities.

### 8.2 Proper motion

The proper motions of objects in SDSS are found automatically by matching the stellar positions to those in the astrometric catalogue USNO-A V2.0 (Monet et al., 1998) derived from the digitized Palomar Observatory Sky Survey (POSS) plates, for objects bright enough to have a POSS counterpart. The normalized distribution of the proper motion of spectroscopically confirmed QSOs is shown in figure 12 as a thick line. These alleged proper motions of QSOs are obviously due to measurement errors, both from POSS and SDSS (primarily from the former). If both measurements follow bivariate gaussian errors with radially symmetric variances $`\sigma _1^2`$ and $`\sigma _2^2`$, then the distribution of the separation of the measurements (i.e., proper motion) also follows a bivariate gaussian with variance $`\sigma ^2=\sigma _1^2+\sigma _2^2`$. Thus, the probability density distribution of observing a star with true proper motion $`\mu _0`$ at an observed proper motion $`\mu `$ is

$$\varphi _{\mu _0}(\mu )=\frac{1}{2\pi \sigma ^2}_0^{2\pi }\mathrm{exp}\left(\frac{\mu ^2+\mu _0^22\mu \mu _0\mathrm{cos}\theta }{2\sigma ^2}\right)\mu 𝑑\theta =\frac{\mu }{\sigma ^2}\mathrm{exp}\left(\frac{\mu ^2+\mu _0^2}{2\sigma ^2}\right)I_0\left(\frac{\mu \mu _0}{\sigma ^2}\right)$$
(4)
where $`I_0(x)`$ is the modified Bessel function of the first kind. Since $`\mu _00`$ for the QSO distribution, the second moment of the distribution of $`\mu `$ should equal $`2\sigma ^2`$ for gaussian statistics. This determination yields $`\sigma _{\mathrm{QSO}}=4.77\mathrm{mas}\mathrm{yr}^1`$, and the resulting theoretical error distribution function from equation 4 is shown in figure 12 as a thin line.

We have estimated the statistical intrinsic proper motions of our sample using the derived photometric distances. The transverse velocity of each star was drawn 1000 times at random from an isotropic velocity distribution with dispersion $`100\mathrm{km}\mathrm{s}^1`$ in one dimension and corrected for an assumed solar velocity $`𝐯_{}=(10,225,7)\mathrm{km}\mathrm{s}^1`$ (Dehnen & Binney, 1998). These choices were based on the kinematic results presented in Paper II. The distribution of intrinsic proper motions so derived, $`N(\mu _0)`$, is shown in figure 12 as the dotted line. The theoretical proper motion distribution for BHB stars with measurement errors should follow the form

$$\varphi _{\mathrm{BHB}}(\mu )_0^{\mathrm{}}N(\mu _0)\varphi _{\mu _0}(\mu )𝑑\mu _0$$
(5)
and is drawn in figure 12 as the dashed line. The distribution of *measured* proper motions of our BHB sample is also shown. This distribution is clearly not very much broader than that of the QSOs; the second moment yields $`\sigma _{\mathrm{BHB}}=4.84\mathrm{mas}\mathrm{yr}^1`$, which is actually too small to account for the combination of measurement error ($`\sigma _{\mathrm{QSO}}=4.77\mathrm{mas}\mathrm{yr}^1`$) and intrinsic proper motion ($`\sigma _{\mathrm{BHB}\mathrm{int}}=2.32\mathrm{mas}\mathrm{yr}^1`$)! In contrast, the proper motion distribution of a sample of faint (here $`r>16.7`$) F stars from SDSS is also shown, with a much broader distribution. This finding demonstrates that contamination in the BHB sample from F stars must be very small. We estimate that the F star contamination is less than $`5\%`$ among the 837 stars in our sample that were matched to USNO-A. Without better understanding of systematic errors, we cannot claim to have detected proper motion of the BHB sample. However, the USNO-B catalog, which has appreciably better data quality than USNO-A, will soon be provided as the first epoch for proper motions, so we defer further analysis to later work.

### 8.3 A map of the Galactic halo and the Sagittarius dwarf tidal stream

Figure 13 shows two views of a three-dimensional representation of the locations of the BHB stars. Each star is colored according to its line-of-sight radial velocity with the solar velocity of $`(10,225,7)\mathrm{km}\mathrm{s}^1`$ (Dehnen & Binney, 1998) added. These two views compose a stereograph which can be viewed by focusing the eyes on a point midway to the paper. The three-dimensional positions of the stars are derived from their angular positions and photometric distances (§7), assuming $`R_0=8\mathrm{kpc}`$.

One could search for substructure in the Galactic halo by a careful examination of figure 13, and indeed a clump appears in the upper left part of the plot (at $`(\mathrm{},b)(350,50)`$). This clump belongs to the leading arm of the Sagittarius dwarf tidal stream (Ivezić et al., 2000; Yanny & Newberg et al., 2000; Ibata et al., 2001a, b). It is interesting to investigate this structure in velocity phase space, as shown in figure 14. The Sagittarius stream approximately follows the great circle on the celestial sphere with pole $`(\alpha ,\delta )=(308,58)`$ or $`(\mathrm{},b)=(94,11)`$, as determined from the distribution of SDSS candidate RR Lyrae stars. Practically the same position is implied by the distribution of 2MASS M giants (Majewski et al., 2003). This great circle intersects the footprint of the SDSS in three places. These three intersections with BHB stars (two in the north and one in the south, with a strip width of 10 degrees) are represented in the three panels in figure 14. One can see a clump in the phase space of the northern subsample representing the Sagittarius stream (left plot); indeed it is the same clump in the upper left of figure 13. The coherent structure in the middle panel (at $`D30\mathrm{kpc}`$) is also associated with the Sgr dwarf tidal stream (trailing arm). This is the same clump that was discovered by Yanny & Newberg et al. (2000) using A-colored stars, and is also detected in the velocity-distance distribution of SDSS candidate RR Lyrae stars (Ivezić et al., 2003). The fact that we see such structure in the BHB sample strongly indicates that our BHB absolute magnitudes and selection criteria are reliable.

## 9 Summary

Using broadband colors and spectroscopic surface-gravity indicators, we have selected $`1170`$ blue horizontal branch stars from the SDSS. For the $`733`$ bright stars ($`g<18`$) the contamination is probably $`10\%`$, and for the $`437`$ faint BHB stars ($`g>18`$) the contamination is probably $`25\%`$. The absolute magnitude of BHB stars was determined from globular cluster color-magnitude diagrams, and an absolute-magnitude-color relation established. This was used to estimate photometric distances, with errors $`0.2`$ mag. The stars also have measured radial velocities with errors $`26\mathrm{km}\mathrm{s}^1`$. Proper motions determined by matching against the digitally scanned Palomar Sky Survey (USNO-A, Monet et al., 1998) also indicate that contamination from the disk is small. We emphasize, however, that no kinematic information was used to select the sample. The data can be downloaded as Table 3 in the electronic version of this paper.

We thank A. Gould, B. Paczyński, R. Lupton, and J. Hennawi for helpful discussions, and M. Strauss for a close reading of a draft of this paper. GRK is grateful for generous research support from Princeton University and from NASA via grants NAG-6734 and NAG5-8083. Funding for the creation and distribution of the SDSS Archive has been provided by the Alfred P. Sloan Foundation, the Participating Institutions, the National Aeronautics and Space Administration, the National Science Foundation, the U.S. Department of Energy, the Japanese Monbukagakusho, and the Max Planck Society. The SDSS Web site is http://www.sdss.org/. The SDSS is managed by the Astrophysical Research Consortium (ARC) for the Participating Institutions. The Participating Institutions are The University of Chicago, Fermilab, the Institute for Advanced Study, the Japan Participation Group, The Johns Hopkins University, Los Alamos National Laboratory, the Max-Planck-Institute for Astronomy (MPIA), the Max-Planck-Institute for Astrophysics (MPA), New Mexico State University, the University of Pittsburgh, Princeton University, the United States Naval Observatory, and the University of Washington.

## Appendix A Remark on contamination and incompleteness

There is no independently confirmed subsample of BHB stars, but we can estimate the accuracy of our spectroscopic selection methods by intercomparing their results.

After applying our color cut, we are left with $`N+C=2338`$ candidates bright enough to be subjected to our spectroscopic selection criteria. Of these, $`N`$ are true BHB stars, and the other $`C`$ are potential contaminants. Let $`\eta _i1`$ be the probability that method $`i`$ fails to identify a true BHB as such, and $`\kappa _i^{}`$ be the probability that it mistakenly accepts an individual member of the non-BHB group. Then the expected number of stars accepted by method $`i`$ is

$$N_i=(1\eta _i)N+\kappa _i^{}C$$
(A1)
Similarly, the numbers of stars selected jointly by methods $`i`$ & $`j`$, or by $`i`$, $`j`$, & $`k`$ , are

$`N_{ij}`$ $`=`$ $`(1\eta _i)(1\eta _j)N+\kappa _i^{}\kappa _j^{}C,`$
$`N_{ijk}`$ $`=`$ $`(1\eta _i)(1\eta _j)(1\eta _k)N+\kappa _i^{}\kappa _j^{}\kappa _k^{}C,`$
respectively, assuming that these methods are truly independent. It is convenient to define $`\kappa _i\kappa _i^{}C/N_i`$ so that the number of contaminants among the $`N_i`$ stars accepted by method $`i`$ is $`\kappa _iN_i`$. Then

$`N_i`$ $`=`$ $`{\displaystyle \frac{(1\eta _i)}{(1\kappa _i)}}N,`$ (A2)
$`N_{ij}`$ $`=`$ $`\left[{\displaystyle \frac{(1\kappa _i)(1\kappa _j)}{N}}+{\displaystyle \frac{\kappa _i\kappa _j}{C}}\right]N_iN_j`$ (A3)
$`N_{ijk}`$ $`=`$ $`\left[{\displaystyle \frac{(1\kappa _i)(1\kappa _j)(1\kappa _k)}{N^2}}+{\displaystyle \frac{\kappa _i\kappa _j\kappa _k}{C^2}}\right]N_iN_jN_k`$ (A4)
As described in §4, we have four methods. Therefore we have six double intersections (A3) and four triple intersections (A4) with which to solve for the six unknown parameters $`(\kappa _1,\mathrm{}\kappa _4,N,C)`$ that they contain. So the problem is overdetermined, allowing us perhaps to test the assumption of independence. Perfect consistency is not to be expected since even for fixed values of the above parameters, the numbers $`\{N_i,N_{ij},N_{ijk}\}`$ are subject to statistical fluctuations.

If the contamination fractions $`\kappa _i`$ are reasonably small—it is sufficient that $`\kappa /(1\kappa )\sqrt{C/N}`$—then the terms involving $`C`$ in eqs. (A3)-(A4) are unimportant, so that

$$1\kappa _i\frac{N_{ij}N_{ik}}{N_{ijk}N_i}$$
(A5)
This approach yields the results shown in Table 5, and also an estimate for the true number of BHBs tested (*i.e*, not including those too faint for spectroscopic selection): $`N772`$. Note that there are three possible choices of $`(j,k)`$ with which to estimate $`\kappa _i`$ from eq. (A5), all of which are shown in the Table to give an idea of their consistency. It can be seen that the contamination and incompleteness fractions for each are indeed small, $`<10\%`$, if these methods are indeed independent as we have assumed. In support of the latter, the differences among the estimates in columns 2-4 are $`\pm 0.01`$ as would be expected by chance since, for example, binomial statistics predict that the number of false negatives for method $`i`$ should fluctuate by $`\pm \sqrt{\eta _i(1\eta _i)N}`$, or about 7 out of 772 if $`\eta _i=0.07`$.