# The GIRAFFE Inner Bulge Survey (GIBS). III. Metallicity distributions and kinematics of 26 Galactic bulge fields

M. Zoccali, S. Vasquez, O. Gonzalez, E. Valenti, A. Rojas-Arriagada, J.

Minniti, M. Rejkuba, D. Minniti, A. Mcwilliam, C. Babusiaux, et al.

► To cite this version:

M. Zoccali, S. Vasquez, O. Gonzalez, E. Valenti, A. Rojas-Arriagada, et al. The GIRAFFE Inner Bulge Survey (GIBS). III. Metallicity distributions and kinematics of 26 Galactic bulge fields. Astronomy and Astrophysics - A&A, 2017, 599, pp.A12. 10.1051/0004-6361/201629805. hal-02177829

HAL Id: hal-02177829

https://hal.science/hal-02177829

Submitted on 17 Sep 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# The GIRAFFE Inner Bulge Survey (GIBS)

## III. Metallicity distributions and kinematics of 26 Galactic bulge fields ⋆,⋆⋆

M. Zoccali<sup>1,2</sup>, S. Vasquez<sup>1,2</sup>, O. A. Gonzalez<sup>3</sup>, E. Valenti<sup>4</sup>, A. Rojas-Arriagada<sup>2,1</sup>, J. Minniti<sup>1,2</sup>, M. Rejkuba<sup>4,5</sup>, D. Minniti<sup>6,2,7</sup>, A. McWilliam<sup>8</sup>, C. Babusiaux<sup>9</sup>, V. Hill<sup>10</sup>, and A. Renzini<sup>11</sup>

1 Instituto de Astrofísica, Pontificia Universidad Católica de Chile, Av. Vicuña Mackenna 4860, 782-0436 Macul, Santiago, Chile
e-mail: mzoccali@astro.puc.cl

2 Millennium Institute of Astrophysics, Av. Vicuña Mackenna 4860, 782-0436 Macul, Santiago, Chile

3 UK Astronomy Technology Centre, Royal Observatory, Blackford Hill, Edinburgh, EH9 3HJ, UK

4 European Southern Observatory, Karl-Schwarzschild Strasse 2, 85748 Garching, Germany

5 Excellence Cluster Origin and Structure of the Universe, Boltzmannstr. 2, 85748 Garching bei München, Germany

6 Departamento de Fisica, Facultad de Ciencias Exactas, Universidad Andres Bello, Av. Fernandez Concha 700, Las Condes, Santiago, Chile

7 Vatican Observatory, V00120 Vatican City State, Italy

8 The Observatories of the Carnegie Institution of Washington, 813 Santa Barbara St., Pasadena, CA 91101-1292, USA

9 GEPI, Observatoire de Paris, CNRS UMR 8111, Université Paris Diderot, 92125 Meudon Cedex, France

10 Université de la Côte d'Azur, Observatoire de la Côte d'Azur, CNRS, Laboratoire Lagrange, Bd de l'Observatoire, CS 34229, 06304 Nice Cedex 4, France

11 INAF-Osservatorio Astronomico di Padova, vicolo dell'Osservatorio 5, 35122 Padova, Italy

Received 28 September 2016 / Accepted 27 October 2016

## ABSTRACT

**Context.** Several recent studies have demonstrated that the Galactic bulge hosts two components with different mean metallicities, and possibly different spatial distribution and kinematics. As a consequence, both the metallicity distribution and the radial velocity of bulge stars vary across different lines of sight.

**Aims.** We present here the metallicity distribution function of red clump stars in 26 fields spread across a wide area of the bulge, with special emphasis on fields close to Galactic plane, at latitudes $b = -2°$ and $b = -1°$, that have not been explored before.

**Methods.** This paper includes new metallicities from a sample of approximately 5000 K giant stars, observed at spectral resolution $R \sim 6500$, in the Calcium II Triplet region. These represent the main dataset from the GIRAFFE Inner Bulge Survey. As part of the same survey we have previously published results for a sample of approximately 600 K giant stars, at latitude $b \sim -4°$, derived from higher resolution spectra ($R = 22\,500$).

**Results.** The combined sample allows us to trace and characterize the metal poor and metal rich bulge populations down to the inner bulge. We present a density map for each of the two components. Contrary to expectations from previous works, we found the metal poor population to be more centrally concentrated than the metal rich one, and with a more axisymmetric spatial distribution. The metal rich population, on the other hand, is arranged in a boxy distribution, consistent with an edge-on bar. By coupling metallicities and radial velocities we show that the metal poor population has a velocity dispersion that varies rather mildly with latitude. On the contrary, the metal rich population has a low velocity dispersion far from the plane ($b = -8.5°$), yet has a steeper gradient with latitude, becoming higher than the metal poor one in the innermost field ($b = -1°$).

**Conclusions.** This work provides new observational constraints on the actual chemodynamical properties of the Galactic bulge, that will help discrimination between different formation models.

**Key words.** stars: abundances – Galaxy: bulge – Galaxy: structure – Galaxy: kinematics and dynamics – Galaxy: stellar content

## 1. Introduction

Until a decade ago, our knowledge about the properties of the Galactic bulge, the second most massive component of the Milky Way, was relatively poor, yet we thought we knew most of what was relevant. We knew it was an elongated spheroid (i.e., a bar,

Stanek et al. 1994), with a broad metallicity distribution (Rich 1988) containing mostly old stars (Ortolani et al. 1995). The results of the investigations cited above were refined by several other authors (e.g., McWilliam & Rich 1994; Zoccali et al. 2003) but none of them covered a wide area of the bulge, thus the common belief was that the bulge was a rather uniform population.

Wide area photometric and spectroscopic surveys conducted in the last ten years revealed that the bulge is much more complex.

Concerning the bulge three-dimensional structure, we now know that the bar flares up in a boxy/peanut, or X-shape

<footnote>* Based on observations taken with ESO telescopes at the La Silla Paranal Observatory under programme IDs 71.B-0617, 385.B-0735 and 187.B-0909.</footnote>

<footnote>** Full Table 1 is only available at the CDS via anonymous ftp to cdsarc.u-strasbg.fr (130.79.128.5) or via http://cdsarc.u-strasbg.fr/viz-bin/qcat?J/A+A/599/A12</footnote>
(McWilliam & Zoccali 2010; Nataf et al. 2010; Saito et al. 2011; Wegg & Gerhard 2013), the unequivocal signature of a formation scenario starting from a disk, whose dynamical instabilities funnel stars, and possibly gas, toward the center, forming a bar that later bends and buckles giving rise to the X-shape (Patsis et al. 2002; Athanassoula 2005). We also know, however, that the oldest and more metal poor stars do not follow the same structure. In fact, RR Lyrae trace a more axisymmetric component compared with the RC (Dékány et al. 2013; Pietrukowicz et al. 2015). Also, metal poor, red clump (RC) stars do not trace the X-shape (Ness et al. 2012; Rojas-Arriagada et al. 2014). Furthermore, there is evidence for an axisymmetric structure within the inner 250 pc of the bulge (Gonzalez et al. 2011a; Gerhard & Martinez-Valpuesta 2012) and for a thin extension of the bar, called the *long* bar (Wegg et al. 2015, and references therein). A recent review of the 3D structure of the Galactic bulge can be found in Zoccali & Valenti (2016).

Another fundamental characteristic of any stellar system is the metallicity distribution function (MDF); indeed, stellar systems may be differentiated by their mean metallicity, or [Fe/H], and by the metallicity dispersion. At any moment, the gas-phase metal content consists of the integral of the nucleosynthetic yields of the preceeding history of star formation, resulting in an increase in the mean metallicity with time. The observed MDF provides constraints on models of chemical evolution. For example, the low mean metallicity of stars in the Galactic halo led Hartwick (1976) to conclude that the halo must have lost its gas through outflows before enrichment to higher metallicity could occur. On the other hand, the lack of metal poor stars in the solar neighborhood, relative to a Closed Box of chemical enrichment, is evidence of gas inflow during chemical evolution (e.g., Pagel 1989). Additionally, the mean MDF may reveal radial or vertical gradients, possibly indicating viscous flows, the radial action of bars, or dissipational collapse.

As mentioned above, the fact that the bulge has a broad MDF has been known for a long time. However, only in the last ~15 years has the measured MDF been precise enough to allow important conclusions to be made on the bulge properties. Zoccali et al. (2008), Hill et al. (2011), Ness et al. (2013a), Rojas-Arriagada et al. (2014), Gonzalez et al. (2015) demonstrated that the bulge MDF is bimodal, with the two peaks having metallicity a few dex below and a few dex above solar, respectively. While the mean metallicity changes across the bulge area, the position of the two peaks does not, thus the mean metallicity gradient, also known for a long time from photometry or low resolution spectroscopy (e.g., Terndrup 1988; Minniti et al. 1995), is in fact due to a different relative fraction of the two metallicity components. None of the studies mentioned above observed stars at latitudes $|b| < 4^\circ$. Due to the higher extinction close to the plane, previous studies of metallicitites in the inner bulge ($|b| < 4^\circ$) are either based on very low number statistics, that is, one or two dozen stars (e.g., Rich et al. 2007, 2012; Schultheis et al. 2015), or on approximately 100 stars with poorer metallicities, compared with those presented here (Babusiaux et al. 2014).

Interestingly, the bulge metal-rich and metal-poor sub-populations have vertical scale heights reminiscent of the local thin and thick disks, respectively, consistent with the $[\alpha/Fe]$ ratios of the thin and thick disks (Alves-Brito et al. 2010; Hill et al. 2011; Gonzalez et al. 2011b). However, these bulge sub-populations have ~0.4 dex higher [Fe/H] than the average thin and thick disk stars in the solar neighborhood. Thus, if these bulge sub-populations truly represent the inner thin and thick disks, a radial [Fe/H] gradient of approximately -0.05 dex/kpc

is implied for both disks (McWilliam 2016). If the bulge was built by the growth and buckling of a stellar bar, then inner thin and thick disk stars entrained into the bulge must have retained vertical scale heights characteristic of their origin, and resulting in the present-day vertical metallicity gradient. The MDF and the chemical composition of bulge stars have recently been reviewed by Ness & Freeman (2016) and McWilliam (2016).

Finally, an observational tool that can be critical in discriminating between different formation scenarios is the stellar kinematics. Indeed, galactic bars usually present streaming motions along the major axis, due to the elongated orbits that sustain the bar. For the same reason, the phase space of stars in the bar shows a significant vertex deviation, contrary to the stars in a spheroidal component, showing a more isotropic velocity distribution with a higher dispersion. Although Saha & Gerhard (2013) demonstrated that a pre-existing spheroid would increase its rotation velocity if a bar is formed at later times, the other kinematical differences mentioned above help discriminate between two or more populations if they belong to different spatial structures. The interested reader may refer to the review on bulge kinematics by Babusiaux (2016). The most relevant of recent investigations demonstrated that the bulge displays cylindrical rotation (Howard et al. 2009; Ness et al. 2013b; Zoccali et al. 2014), with a velocity dispersion increasing toward the Galactic plane, including a peak in the inner ~2 degrees from the Galactic center (Zoccali et al. 2014, hereafter Paper I). The kinematical properties of the metal poor and metal rich bulge components are slightly different, with the metal poor having higher velocity dispersion, at least in the outer bulge (Ness et al. 2013b).

We present here the MDF derived for 26 bulge fields, observed within the GIRAFFE Inner Bulge Survey (GIBS), an ESO Large Programme (ID 187.B-0909; PI: Zoccali) carried out with the GIRAFFE spectrograph of the FLAMES instrument (Pasquini et al. 2002) at the Very Large Telescope (VLT). The aim of the survey is to investigate how the metallicity and radial velocity distributions vary across the bulge area, with close attention paid to the fields close to the Galactic plane, at $b = -2$ and $b = -1$, not explored by other surveys. Here we combine the metallicities of the stars with their radial velocities, already discussed in Paper I, in order to characterize the spatial distribution and kinematics of the two bulge metallicity components.

## 2. The data

The present investigation is based on spectra for approximately 5500 RC stars in several fields across the Galactic bulge, as shown in Fig. 1. All of them are located within the region mapped by the VVV Survey ($-10 < l < +10$ and $-10 < b < +5$; Minniti et al. 2010).

The targets were selected from the near infrared CMDs obtained within the VVV Survey, and calibrated to the 2MASS photometric system as in Gonzalez et al. (2011c). The target selection, discussed in Paper I, is shown in Fig. 2: the GIRAFFE targets are RC stars, in a relatively narrow range of apparent magnitude around the RC peak. The size of the target box is allowed to vary in different fields, according to the stellar density of the field, and the number of targets observed, spanning on average $\Delta J_0 \approx 0.5$ mag, corresponding to approximately 1 kpc along the line of sight. In Fig. 2 the observed stars are color-coded according to their metallicity (see caption), and a color histogram of the targets is compared with that of all the RC stars. This is meant to emphasize that our target selection does not impose any significant bias in metallicity, except excluding the most metal poor stars ($[\text{Fe/H}] \le -1.5$) that, if older than 10 Gyr,
<figure><img src="image_2.png" /></figure>

**Fig. 1.** Position in the sky of the 26 GIBS fields presented here. The large square is the area covered by the VVV Survey. Red circle are fields observed at low resolution (LR8), purple squares are fields observed at high resolution (HR13; Paper II). The blue star is the Baade's Window field used for the CaT calibration (Vásquez et al. 2015), and the blue triangle is the (0, -6) field discussed in Zoccali et al. (2008).

would end up in the blue horizontal branch. The absence of a metallicity bias is especially relevant in the interpretation of the MDFs for the two innermost fields (one of which shown here) where we find a significant metal poor population.

Two independent sets of stars were observed within each field, for a total of approximately 210 stars. The two sets are identified as F1 and F2 in our tables. In a few fields that we considered particulary interesting, four sets of stars were observed instead, for a total of approximately 430 stars. Those fields are the two closest to the Galactic center, at (0, -1) and (0, -2), and the field at (0, -8) because it shows the double RC feature, signature of the presence of X-shape (or boxy/peanut shape) as demonstrated by McWilliam & Zoccali (2010), Nataf et al. (2010), Saito et al. (2011), Wegg & Gerhard (2013). In the (0, -8) field, we observed two sets of stars in the bright RC (bRC-F1 and bRC-F2) and two sets in the faint RC (fRC-F1 and fRC-F2).

A detailed description of the target selection strategy, observations and spectrum extraction and calibration is given in Paper I. Table 1 in that paper lists the coordinates of all the fields, together with the number of targets, the spectral setup, and the total exposure time associated to each of them. Here we recall that most of the GIBS fields were observed with the LR8 setup, covering the Calcium II triplet (CaT) feature (8206–9400 Å) at resolution R = 6500. Four fields at b ~ -3.5° were observed at higher resolution using setup HR13 (6120–6405 Å). Abundances for those are discussed in Gonzalez et al. (2015, hereafter Paper II).

Spectra for stars in two extra fields were added here to the original Large Programme. The first one consists of a sample of 178 RC stars in Baade's Window, at (l, b) = (1, -4), marked as a star in Fig. 1. These stars, observed through setup HR13, were presented in Hill et al. (2011) but have been re-analyzed in Paper II, and we adopted here metallicities from the latter work. A subset (111) of these stars were also observed with setup LR8, and 80 of them were used to derive the CaT to [Fe/H] calibration, presented in Vásquez et al. (2015) and adopted here. The second additional dataset consists of 213 red giant branch (RGB) stars in a field at (l, b) = (0, -6) from Zoccali et al. (2008). Note that these stars are not the same as

<figure><img src="image_3.png" /></figure>

**Fig. 2.** Example of the GIBS target selection criteria for the LRp0m1 field. Top: the VVV CMD together with the 432 GIRAFFE targets, color-coded according to their metallicity (blue: [Fe/H] < -0.5 dex: green: -0.5 < [Fe/H] < -0.1; orange: -0.1 < [Fe/H] < +0.3 and red: [Fe/H] > +0.3.). Bottom: color histogram of all the stars in the FLAMES field, in the black box of the CMD, compared with the color histogram of the selected targets.

included in Paper I. Because it was focused on kinematics only, Paper I included a sample of 454 RC stars, at (0, -6), from Vásquez et al. (2013), that are more numerous than the sample analyzed here. However, those spectra had lower S/N and, most importantly, the dispersion direction encompassed the four detectors of the IMACS@Magellan mosaic, preventing a reliable continuum definition and therefore hampering the metallicity determination. With the focus on the MDF, here we replaced the Vásquez dataset with the RGB stars from Zoccali et al. (2008). The latter, although smaller in number, have been observed with the same instrument as our HR sample. The analysis carried out in Zoccali et al. (2008) is compatible with the one in Paper II within 0.1 dex, as verified in the Baade's Window field. The targets at (l, b) = (0, -6) are the only ones in the present work that lie on the RGB, slightly above the RC.

In summary, we discuss here the MDF for 26 fields, shown in Fig. 1. For 20 of them, shown as red circles, we have LR8 spectra with R = 6500. Metallicities were derived by measuring CaT equivalent widths and converting them to [Fe/H] following the recipe in Vásquez et al. (2015). For six other fields (blue squares, star and triangle) we have spectra at R = 22500 and [Fe/H] measurements were derived from equivalent widths of isolated FeI lines, constraining also stellar surface parameters.

## 3. Iron abundances

### 3.1. *High-resolution spectra*

The analysis of the high-resolution spectra was presented in Paper II, where the details of the adopted method to obtain [Fe/H] and [Mg/Fe] can be found. Here we briefly recall how [Fe/H] abundances were derived.

Stellar atmospheric parameters were discussed in Paper II; namely effective temperature (T<sub>eff</sub>), surface gravity (log g), microturbulent velocity (v<sub>μ</sub>), and metallicity ([Fe/H]) using the
<table><caption>**Table 1.** Coordinates, VVV magnitudes, radial velocities, and metallicities of the program stars.</caption><thead><tr><th>ID</th><th>RA (J2000)</th><th>Dec (J2000)</th><th>J</th><th>K<sub>s</sub></th><th>RV</th><th>Err RV</th><th>[Fe/H]</th></tr></thead><tbody><tr><td>LRp0m1_F1_OGLE5_140424</td><td>17:50:20.66</td><td>-29:58:42.30</td><td>14.617</td><td>13.234</td><td>-42.9</td><td>2.1</td><td>-0.33</td></tr><tr><td>LRp0m1_F1_OGLE5_141397</td><td>17:50:19.03</td><td>-29:58:04.10</td><td>14.761</td><td>13.361</td><td>198.0</td><td>2.3</td><td>-0.65</td></tr><tr><td>LRp0m1_F1_OGLE5_145942</td><td>17:50:19.62</td><td>-29:56:20.40</td><td>14.166</td><td>12.551</td><td>253.0</td><td>2.4</td><td>-0.10</td></tr><tr><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td></tr></tbody></table>

**Notes.** The full table is available at the CDS.

standard iterative method based on the equivalent widths of isolated FeI lines. Initial abundances are derived from each line by adopting first-guess photometric temperature and gravity, and these are later refined by imposing excitation and ionization equilibrium. The code GALA (Mucciarelli et al. 2013) has been used here to perform this task automatically.

First-guess photometric temperatures were calculated by applying the calibration by Ramírez & Meléndez (2005) to the ($V-K_s$) magnitudes from OGLEII (Udalski et al. 2002) and the VVV catalogues in Gonzalez et al. (2012). The reddening maps in the latter paper were also used to correct for extinction. Absolute V-band magnitudes were derived using mean line-of-sight distances to each specific field from Gonzalez et al. (2013) and the bolometric corrections from Alonso et al. (1999). These are then used to estimate photometric gravities based on the classical formula:

$$ \log(g) = \log(g_{\odot}) + \log\left(\frac{M_{*}}{M_{\odot}}\right) + 0.4(M_{\text{Bol},*} - M_{\text{Bol},\odot}) + 4\log\left(\frac{T_{\text{eff},*}}{T_{\text{eff},\odot}}\right), \quad (1) $$

where $M_{\text{Bol},\odot} = 4.72$, $T_{\text{eff},\odot} = 5770$ K, and $\log(g_{\odot}) = 4.44$ dex. As usual, we adopted a fixed value of $M_{*} = 0.8 M_{\odot}$. Global metallicity and microturbulent velocity were fixed to 0.0 and 1.5, respectively, as a first step.

An initial ATLAS9 stellar model atmosphere was constructed using these first-guess values. This was later updated by GALA while iteratively searching for spectroscopic effective temperatures and microturbulent velocity, by imposing excitation equilibrium and the null slope of iron abundance versus equivalent width of the Fe line, respectively. As discussed in Gonzalez et al. (2015), the mean uncertainty in the metallicity measurements is $\sigma_{[\text{Fe/H}]} = 0.2$ dex, reaching as low as 0.1 dex at the low metallicity end, and up to 0.4 dex at the high metallicity end.

The method described above has been used identically to derive the abundances for the four GIBS fields observed at high resolution (blue squares in Fig. 1) and for the Baade's Window giants used by Vásquez et al. (2015) to derive the CaT to [Fe/H] calibration. The latter is applied here to the low resolution spectra, thus insuring that the final metallicities adopted for our stars are all in the same scale.

## 3.2. Low resolution spectra

The pre-reduction of the LR8 spectra, together with the fitting of the three CaT lines is described extensively in Vásquez et al. (2015) where the CaT versus [Fe/H] calibration was derived. We summarize the main steps here. The spectra were de-biased, flat-fielded extracted and wave calibrated using the GIRAFFE pipeline (version 2.8) provided by ESO. Sky subtraction was instead applied with IRAF, by first combining the approximately 20 sky fibres in a master sky, which was then subtracted from

<figure><img src="image_4.png" /><figcaption>**Fig. 3.** Top panels: trend of the EW of line C versus line B for two different fiber allocations in the field at (0, -8). The orange points are 103 targets in the bright RC, while the purple ones are 102 targets in the faint RC. The green line is the same trend for the RC stars in Baade's Window used to derive the CaT vs. [Fe/H] calibration. While the orange targets beautifully follow the expected ratio, the purple points have an offset, that might be due either to an overestimation of line B and/or to an underestimation of line C. Bottom panels: (left) MDF derived from the 103 targets with good line ratio; (right) the MDF derived from the targets with offset line ratio, in purple, is compared with the orange one, obtained from targets with good line ratio.</figcaption></figure>

the science spectra by means of the task skytweak. The spectra were then normalized with the IRAF task continuum, and finally cross-correlation (task fxcor) was performed against a template synthetic spectrum in order to bring the spectra to the restframe velocity.

Equivalent widths (EWs) of the two strongest CaT lines (8542 Å and 8662 Å) were measured by fitting the observed flux with the sum of a Gaussian component, fitting the line core, plus a Lorentzian function fitting the wings. Vásquez et al. (2015) presented several tests of the stability of such fits against metallicity and signal to noise (S/N), by means of synthetic spectra. They demonstrated that the weakest CaT line at 8498 Å (hereafter line A), does not behave smoothly with metallicity, as is the case for the two stronger lines (B and C), and it is the most affected by sky subtraction, due to the presence of two sky lines close to its wings. For this reason the CaT-to-[Fe/H] calibration presented in Vásquez et al. (2015) relies only on the two stronger lines.

Figure 5 in Vásquez et al. (2015) demonstrates that, for the stars in Baade's Window, used to derive the CaT-to-[Fe/H] calibration, the ratio between line C and line B is remarkably constant along metallicity. The ratio of these two lines was compared with that of Baade's Window for the stars in all the fields analyzed here. Most of the stars indeed follow the same line ratio
<figure><img src="image_5.png" /><figcaption>Fig. 4. Metallicity distribution function (MDF) for the 26 fields of the present study. Orange histograms have been derived from CaT spectra, while green histograms come from the high-resolution spectra discussed in Paper II. Blue histograms are from our previous works. Panels are arranged in the same way as the fields in the sky, with approximate coordinates indicated in the labels. The number of stars observed in each field is also labeled for each field.</figcaption></figure>

as the calibrators in Baade’s Window. In some fields, however, all the stars observed within a given exposure show a shift in the ratio of line C to line B, with respect to Baade’s Window (Fig. 3). The shift is constant for all the stars, that is, it does not depend on the metallicity nor the surface parameter of the individual stars. It is also independent from the radial velocity of each star, which differs by up to 400 km s⁻<fnref n="1" />. We conclude that it must be due to an improper continuum normalization (i.e., a residual slope in the pseudo-continuum of the two lines). A variation of less than 1% in the EWs (upper right panel of Fig. 3) corresponds to a shift in the pseudo-continuum that is very hard to detect. Despite our many attempts, it was impossible to relate the presence of the offset to any specific sky condition (Moon, humidity, etc.), nor to correct it. Fortunately, the impact of the offset on the final metallicity is very small, if any (lower right panel). We therefore did not correct for this effect, keeping in mind that it happens, with the same magnitude shown in Fig. 3, for 10 sets of stars, in 7 fields out of the 20 observed with the CaT setup<fnref n="1" />.

In order to apply the CaT-to-[Fe/H] calibration equation from Vásquez et al. (2015) one needs to calculate the so-called reduced EW parameter, W', defined as

$$W' = \Sigma EW + 0.384(K - K_{RC}), \quad (2)$$

where ΣEW is the sum of the EWs of line B and line C, and the term 0.384(K − KRC), is a correction that simultaneously takes into account the effect of the surface gravity and temperature of the star, parametrized as the magnitude difference between the star and the RC in that specific field.

To that end, the observed Ks magnitudes of all stars were corrected using the interstellar extinction maps of Gonzalez et al. (2012) before calculating the difference with respect to peak reddened RC magnitude in that field. Since the stars were chosen to be on the RC, this difference is very close to zero for all our targets. We nevertheless included this correction because the target selection was not carried out on the observed magnitudes before the extinction map was published. Therefore, in the innermost, most extincted fields, some stars could end up slightly displaced from the peak RC in the de-reddened CMD. The CaT calibration presented in Vásquez et al. (2015) has a rms uncertainty of 0.2 dex.

The coordinates, magnitudes, metallicities and radial velocities of the GIBS targets discussed here are given in Table 1.

## 4. The bulge Metallicity Distribution Function

The resulting metallicity (iron) distribution function for the 26 fields is shown in Fig. 4. Panels in this figure are arranged in the same way as the fields in Fig. 1, with labels within each panel showing the galactic coordinates of each field, and the number

<footnote n="1">¹ The affected sets of stars are the following: LRp5p4-F1, LRp5p4-F2, LRp0m1-F3, LRp3m2-F2, LRp0m2-F2, LRp0m2-F3, LRm8m2-F2, LRp0m8-fRC-F2, LRm3m8-F2, LRm8m8-F2.</footnote>
of stars included in the histogram. The orange histograms refer to the fields observed at low resolution (Sect. 3.2), the green histograms are the fields observed at high resolution (Sect. 3.1) with our ESO Large Programme, and the blue histograms are the 213 RGB stars from Zoccali et al. (2008) at (0, −6), and the 178 RC stars from Hill et al. (2011) in Baade’s Window.

All the histograms presented in this Section have been obtained as Average Shifted Histograms. Specifically, three histograms were constructed, each with a bin of 0.21 dex and starting points at [Fe/H] = −2.0, −1.93 and −1.86 (i.e., shifted by 0.07 dex), respectively. A vertical line, fixed at [Fe/H] = 0, has been drawn in each panel as a guide to the eye. All the fields in Fig. 4 show a clear bimodality, apart from two of them, LRp8m2 and LRm8m6, that appear unimodal. In order to assess the significance of the bimodality, we ran a Gaussian Mixture Model (GMM Muratov & Gnedin 2010) allowing for the presence of one, two, or three Gaussians, with different means and dispersions. The best fitting solution for most of the fields is the trivariate one, with a very small peak for [Fe/H] < −0.8. This peak usually constitutes a handful of stars, and in some fields is so small that it is not even visible in Fig. 4. Suppressing this peak, and fitting only two Gaussians, however, resulted in the metal-poor Gaussian being very broad; enough to include these few stars, and therefore extending to supersolar metallicities. This peak would unrealistically include some of the most metal rich stars of the sample. The triple Gaussian fit did not converge in four fields (LRp5p4, LRm5m2, HRp8m3 and Baade’s Window); for these fields the data were best fit with double Gaussians.

In fields LRp8m2 and LRm8m6, the two-Gaussian model did not yield a significantly better fit than a single Gaussian. It would seem unrealistic that the bulge contains two metallicity components, with a separation of ~0.6 dex in [Fe/H], except for in these two fields. We therefore double-checked the target selection and the velocity distributions of stars in these two fields, to ensure that we are sampling the same populations as in the others. From Fig. 10 below, we see that at (+8, −2) and (−8, −6), the velocity versus metallicity scatter plots do not show any strange feature, that could indicate a possible sub-population.

Figure 5 (left panels) shows the *p*-value of Kolmogorov-Smirnov test comparing the MDF of these fields to the global MDF for the other fields at the same latitude. In the case of the (+8, −2) field, the *p*-value (9%) does not permit rejection of the hypothesis that this field is drawn from the same MDF. The top right panel shows the target selection for this field: while the selection was done on the observed CMD (before the VVV reddening maps were available), de-reddened magnitudes are shown here; therefore, the targets are no longer confined to a CMD box. This field in particular is very close to the plane, where differential reddening is non negligible. After correction, some targets seem to belong to the foreground disk rather than to the bulge. Targets here are color coded according to their metallicity: blue points are stars with [Fe/H] < −0.6 dex, green with −0.6 < [Fe/H] < −0.1 dex, orange with −0.1 < [Fe/H] < +0.3 dex, and red with [Fe/H] > +0.3 dex. Contaminating disk stars all belong to the peak at [Fe/H] ~ −0.2 dex. Indeed, removing stars with color (J − K)₀ < −0.35 would change the *p*-value to 15%. We conclude that this field happened to have a larger disk contamination due to the higher differential reddening. Its distribution is nonetheless compatible, within the statistical fluctuations of the present sample, with the global one at this latitude. Note also that the metallicity derived for disk stars would not be accurate, since their gravity differs significantly from that of CaT calibrating stars.

<figure><img src="image_6.png" /><figcaption>**Fig. 5.** Left panels: MDF of the discrepant fields, compared with the global MDF of the other fields at the same latitude. The *p*-value of the Kolmogorov-Smirnov test comparing these peculiar MDFs to the global MDF is indicated. Right panels: corresponding target selection for these fields. Targets are color-coded according to their metallicity, with blue dots having [Fe/H] < −0.5 dex, green −0.5 < [Fe/H] < −0.1, orange −0.1 < [Fe/H] < +0.3 and red [Fe/H] > +0.3.</figcaption></figure>

Stars in the field at (−8, −6) are shown in the middle panels. Here the MDF is much more metal poor than the global one, with null probability (10⁻¹¹) of having been drawn from it. On the other hand, the CMD does not show obvious contaminations, and neither does the velocity versus metallicity plot in Fig. 10. We do not have a plausible explanation for the observed difference between the MDF of this field and the global one at this latitude. Even if a unimodal Gaussian distribution would fit the MDFs of these two fields marginally better than a bimodal one, a bimodal fit was imposed for consistency with the other fields. Obviously the position of the peaks will not be well constrained in this case.

Finally, another rather peculiar field is the one at (+5, +4). In this case the MDF has two main peaks, similarly to all the others, however the entire range is shifted to metal rich values. Again, the target selection does not seem to differ from that of the other fields. This is one of the fields where the ratio of line B to line C has an offset with respect to the calibrating stars in Baade’s Window. We have shown in Fig. 3 that the impact of this offset on the MDF is negligible, but in this field we have an offset in both target datasets (LRp5p4-F1 and LRp5p4-F2), therefore the impact of the offset on the MDF cannot be properly verified. We consider this field similarly to the others, keeping in mind its peculiarity and ensuring that none of the conclusions drawn in the present paper are based on any of the fields in Fig. 5. Rather, we draw robust results in spite of the peculiar shape of the MDF in those fields.

In the following we attempt to characterize the two metallicity populations seen in the bulge. We assume that the very small peak at [Fe/H] ~ −1 is made up by contaminating halo stars. Notice that our MDFs do not include very metal poor stars because our target selection criterion excludes stars outside the
<figure><img src="image_7.png" /><figcaption>**Fig. 6.** Peak metallicity of the metal-poor (blue) and metal-rich (red) Gaussians fitted in Fig. 4. The names of the fields are listed on the left, from the closest to the plane at the top, to the ones in the outer bulge at the bottom.</figcaption></figure>

color range of the RC. Stars with [Fe/H] < -1, if they are old enough, would become the blue of the RC, and would therefore fall outside the target selection box. The MDFs for the fields at |b| > 3° confirm previous findings by our group (Zoccali et al. 2008), the ARGOS survey (Ness et al. 2013a) and the Gaia-ESO survey (Rojas-Arriagada et al. 2014). Namely, that the bulge has two metallicity components: with the metal-poor one dominating at high latitude, and the metal rich dominating closer to the plane. The mean metallicity, therefore, has a gradient that is due to the different mix of the two populations, and it is mostly vertical rather than radial (Paper II). Figure 6 shows the metallicity of the metal poor (blue) and metal rich (red) peaks as fitted above. Error bars were calculated by means of a non-parametric bootstrap (drawing from the input sample with repetitions) of 100 realizations (Muratov & Gnedin 2010). Overall, the metal poor and metal rich peaks stay constant within the errors, with only two fields showing relative discrepancy, both of them discussed above.

The most important result from the present data is that, at latitude smaller than |b| = 3°, the metal poor component becomes important again. Because the variation of the MDF with longitude, at constant latitude, is smaller than the variation as a function of height from the plane, we have co-added the stars in all the fields at constant latitude, and show the results in Fig. 7. Thanks to the larger statistics, the MDFs here are much smoother and the two Gaussian fits to stars with [Fe/H] > -0.8, obtained with the GMM, are excellent. The label in each panel states the fraction of metal poor stars, that is, the relative contribution of the blue Gaussian to the total bimodal fit. This fraction is 0.77 at b = -8.5°, drops continuously when moving closer to the plane, reaching a minimum of 0.29 at b = -3.5°, and then rises back up to 0.6 at b = -2° and 0.53 at b = -1°.

The MDF for the innermost ~2 degrees from the Galactic center was derived by Rich et al. (2007) and Rich et al. (2012). They analysed near-IR, high-resolution spectra for ~15 M giants in each of two fields at (0, -1) and (0, -1.75), respectively, and found only metal rich stars, with only a few of them having

<figure><img src="image_8.png" /><figcaption>**Fig. 7.** MDFs of all the stars at constant latitudes, where fields at different longitudes have been combined. The fraction of metal poor stars compared to the total is listed in each panel. As discussed in the text, the fraction of metal poor stars decreases from $b = -8.5$ to $b = -3.5$, where it reaches a minimum, but then rises again at latitude closer to the Galactic plane. The vertical dashed lines mark the limits of the metal poor and metal rich populations (see text).</figcaption></figure>

[Fe/H] < -0.3, and none with [Fe/H] < -0.5. Babusiaux et al. (2014), on the other hand, derived the MDF from CaT spectra for 107 giants in a field at (0, +1). They do find a metal poor population, rather similar to the present one at (0, -1), although with a longer metal poor tail, possibly due to the slightly different target selection box. Finally, the presence of a metal poor component within a radius of 1 degree around the Galactic center was detected by Schultheis et al. (2015) from a sample of 33 M giants observed within the APOGEE data.

Two things are very clear from Figs. 4 and 7: i) the bulge does have two metallicity components; and ii) the spatial distribution of the components is different. In order to investigate the density distribution of the two components we need to convolve the information given in Fig. 4 with what we know about the bulge global density distribution. In other words, although the metal poor component is dominant in the outer bulge (at $b = -8.5°$) it might not be very conspicuous, simply because the stellar density at these latitudes is rather low. We investigate this specific point in the next section.

**5. Characterizing the two metallicity populations in the bulge**

In this section we will attempt to characterize the spatial distribution of the metal poor and metal rich stars, separately. We will also investigate the radial velocity distribution of the two populations. Let us stress, that although it is clear that these two populations have different properties, this does not necessarily imply that they have different origins. Certainly, the difference
in spatial and kinematical properties between the two metallicity
components must be taken into account when studying bulge for-
mation models that should ultimately provide information about
the physical processes behind them. The purpose of this paper,
however, is to present and describe the observational data only.

5.1. *Density distribution of red clump stars*

Using the fact that the bulge MDF can be parametrized by two
Gaussian components, each having a coherent spatial variation,
the fraction of metal poor to metal rich stars given by the MDF
in each field can be translated to a stellar density of each com-
ponent individually. This can be done by using the bulge stellar
density map from Valenti et al. (2016), who counted the num-
ber of RC stars across the bulge using PSF photometry from the
VVV survey. We use this number and the relative fraction of
metal-poor and metal-rich Gaussian components of each of the
GIBS fields to calculate the corresponding number of metal poor
and metal rich RC stars, in each field. The number of RC stars
for each component of all fields is listed in Table 3.

In order to visualize the differences in the spatial distribu-
tion of metal poor and metal rich RC stars, Table 3 can be used
to construct spatial density maps for the two components. This
is challenging with the current data because the GIBS field grid
is relatively coarse. In order to improve the spatial resolution,
we folded the field grid both in Galactic longitude and lati-
tude. When doing so, we are imposing that the density maps
must be symmetric both about *l* and *b*. While symmetry about
the Galactic plane is very plausible, the orientation angle of
the Galactic bar does produce an asymmetry in the projected
density between positive and negative longitudes (Alard 2001;
Valenti et al. 2016). For the present investigation, however, we
are assuming that our data are not dense enough to allow us to
detect differences between the near and far sides of the bar.

Figure 8 shows the number of RC stars in two strips at ap-
proximately constant latitude, as a function of longitude, for the
metal poor (top panels, in blue) and the metal rich (bottom pan-
els, in red) components. Filled symbols show the actual positions
of the GIBS fields, at negative (squares) and positive (triangles)
latitude. Open symbols are the mirror positions of these fields,
when folded in longitude. The two discrepant fields LRp8m2 and
LRm8m6 have been omitted from this analysis. Figure 8 shows
that the metal poor component is peaked toward the center, while
the metal rich one has a relatively constant density distribution
at least at these latitudes. The relatively flat distribution of metal
rich stars resembles the RC surface density in Nataf et al. (2015,
their Fig. 13), at |b| = 5.5, although the latter study included all
the RC stars, without distinction on metallicity. Figure 8 shows a
sub-sample of the real data, without any modeling nor fitting. It
has, however, one important limitation, because the fields shown
here have only approximately constant latitude: the fields at |b| ~
4° actually span a range from |b| = 3.5° to |b| = 4.5°. Other fields
not shown here are even more scattered at different coordinates.
A deeper analysis of what the data show, therefore, requires an
interpolation and a two-dimensional visualization.

We applied a bivariate linear interpolation to the irregular
grid of GIBS fields using the *akima* package in R (Akima 1978)
to produce stellar density maps of the metal poor and the metal
rich components, as shown in Fig. 9. To construct the maps, we
normalized the number of RC stars of each component by the
maximum number of metal poor RC stars so that the color code
of the maps can be compared directly. We focus on the general
shape of the maps instead of the details, as the latter are highly
dependent on the interpolation method and the non-uniform grid.

<figure><img src="image_9.png" /><figcaption>Fig. 8. Stellar density of RC stars, assumed to trace the global stellar density, for the metal poor (top panels, blue) and the metal rich (bottom panels, red) components. Filled squares/triangles refer to actual fields at negative/positive latitude, respectively. Open symbols are their mirror positions, when folded in longitude. Two strips at almost constant latitude are shown here, as indicated in the labels. Notice that the metal poor component has a steeper gradient toward the center, while the metal rich one has, if any, a very mild gradient at these latitudes.</figcaption></figure>

For this reason, we smoothed the map using an isotropic Gaus-
sian kernel with standard deviation of 1 deg. The fields LRp8m2
and LRm8m6 have been omitted from the maps because the frac-
tion of metal poor to metal rich components obtained from them
are not reliable in our opinion. In any case, thanks to the four-
fold symmetry that we chose to impose, the maps obtained with
these fields included are only barely distinguishable from the
ones presented here.

It is clear from Fig. 9 that the density map of metal poor
RC stars is more concentrated than that of the metal rich
RC stars. Furthermore, the metal rich map appears more elon-
gated, or “boxy”, along Galactic longitude because of a steep
drop in the stellar counts after |b| = 5 at all Galactic longitudes.
On the other hand, toward the central region of the maps, the
metal poor stars are not only seen to be more centrally concen-
trated but they also overcome the number of metal rich RC stars.

5.2. Kinematics versus metallicity trends

We present here a kinematical characterization of the sample
stars, with particular emphasis on possible differences between
the metal poor and metal rich components. Only the radial ve-
locities, and velocity dispersions, available from the GIBS spec-
tra are discussed here. The proper motions for the spectroscopic
targets are currently under analysis from VVV data. These will
allow us to calculate three-dimensional space velocities, which
will be discussed in a forthcoming paper.

Figure 10 shows the galactocentric radial velocity of the sam-
ple stars versus their [Fe/H], for all the fields. Orange points are
the fields observed with the LR8 setup in the CaT spectral re-
gion. Light blue points are the targets observed at higher spectral
resolution. Green and pink points in the field at (0, −8.5), refer
to stars in the bright and faint RC, respectively, which were also
<figure><img src="image_10.png" /><figcaption>**Fig. 9.** Density map of metal poor (top) and metal rich (bottom) RC stars obtained using the two component MDF in each field and the total number of RC stars from Valenti et al. (2016). Open circles show the grid of fields used in the interpolation, resulting from folding the original map of GIBS fields with respect to both Galactic latitude and longitude. The number of RC stars of each component was normalized by the maximum number of RC stars in the metal poor map so they are color-coded to the same scale between values of 0 and 1. Contours are over-plotted in both maps and correspond to differences of 0.1 in the normalized RC star counts.</figcaption></figure>

observed with the LR8 setup. This figure shows that the distributions are fairly homogeneous. Although there are a few outliers with respect to the main radial velocity distribution, there are no obvious bimodalities, nor tight clumps suggestive of possible streams. Due to the bulge rotation, stars in the fields at positive longitudes (panels on the left side of this figure) have positive (receding) mean velocities, while stars in fields at negative longitudes have negative (approaching) mean velocities.

This figure does not show any obvious trend of the mean radial velocity versus metallicity. It does show, however, that the velocity dispersion has a trend as a function of metallicity.

Comparing, for example, the innermost field (second row from the top) with the outer one along the minor axis (pink and green points), the velocity dispersion of the metal poor stars, smaller than that of the metal rich stars in the innermost field, becomes larger in the outer field. This change is more evident along the minor axis, but is roughly conserved at $l \neq 0$. The black broken lines within each panel show the $\pm 1\sigma$ interval around the median, as a function of metallicity, for three 0.5 dex bins centered at [Fe/H] = -0.7, -0.2, and +0.3, respectively, with the aim of emphasizing the change of behavior. As expected, the stars observed at high resolution at $b \sim -3.5$ (light blue points) show a very noisy behavior when binned in metallicity because they are fewer in number, and also because at this latitude the metal poor component is very small compared to the metal rich one.

The variation of the velocity dispersion discussed above is best illustrated in Figs. 11 and 12, showing the velocity dispersion as a function of longitude, for different latitude, splitting the sample into a metal poor and a metal rich one. With the aim of separating the two components, we have decided to make a constant cut at [Fe/H] = 0.0, with the metal poor sample composed of stars with [Fe/H] < -0.1 dex, while the metal rich only containing stars with [Fe/H] > +0.1 dex. The narrow strip at -0.1 < [Fe/H] < +0.1 has been left out as "no man's land", because the large cross contamination would add more noise than number statistics. We have shown in Fig. 6 that the peaks do not move significantly from one field to the other, hence the choice of a constant cut (e.g., Fig. 7). The GMM yields the membership probability of each star to either population. However, in the overlap region, no model can possibly establish the membership probability on a one by one star basis. Therefore, we rather exclude the region where the overlap between the metallicity distribution is significant.

Figure 11 confirms what was already shown in the ARGOS data (Ness et al. 2013a), that is, that the sigma of metal poor stars is, on average, higher than that of metal rich stars, at least for $|b| > 4^\circ$. It should be noted that the agreement with previous results should be taken qualitatively, because the metallicity cut in different studies has been chosen differently. In addition, the target selection in the ARGOS survey is such that many more metal poor (halo) stars are included. GIBS data allows us to extend the investigation to inner regions, showing that the trend of sigma versus metallicity is inverted at $b = -1^\circ$, at least along the bulge minor axis. Figure 12 best illustrates this point, showing that for $|b| > 3^\circ$ the metal poor component has a higher sigma compared to the metal rich one, with the difference becoming larger moving outward. At latitude $b = -2^\circ$ the two metallicity components have comparable sigma, while at $b = -1^\circ$ it is the metal rich component that shows a higher sigma. Unfortunately we only have one field at this latitude, but there are 432 targets in it, and therefore the velocity dispersion measurement is robust. An indication of the inversion in the sigma versus metallicity trend, at very low latitudes, was already found in the data by Babusiaux et al. (2014 – see also Babusiaux 2016). With its very homogeneous target selection, the GIBS survey allows us to confirm that the inversion is real. This observational result seems to be due to the fact that the dispersion of the metal poor stars increases from the outer bulge ($b = 8^\circ$) toward the plane, although slower than that of metal rich stars, but then stays constant in the three strips at $b = -3.5^\circ, -2^\circ, -1^\circ$. In other words, the velocity dispersion of the metal rich component has a gradient with latitude much stronger than that of the metal poor component.

Kunder et al. (2016) recently suggested that bulge RR Lyrae, tracing the oldest, more metal-poor component, rotate slower than RGB and RC stars. We investigate whether this signature
<figure><img src="image_11.png" /><figcaption>**Fig. 10.** Metallicity versus galactocentric radial velocity for all the targets in the GIBS fields. Orange points are stars observed in CaT,azure points are stars observed at high-resolution. The field at (0, -8) has a double RC: bright-RC targets are shown in green while faint-RC ones are shown in pink. The solid lines mark the ±1σ interval around the mean velocity at three different metallicities. See text for details.</figcaption></figure>

<table><caption>**Table 2.** Slope of $V_{GC}$ vs. longitude, in five latitude strips.</caption><thead><tr><th>Latitude</th><th>[Fe/H] &lt; -0.3</th><th>[Fe/H] &lt; -0.1</th><th>[Fe/H] &gt; +0.1</th><th>[Fe/H] &gt; +0.2</th></tr></thead><tbody><tr><td>b ~ -2°</td><td>8.2 ± 1.4</td><td>9.3 ± 1.1</td><td>10.9 ± 1.0</td><td>10.8 ± 1.1</td></tr><tr><td>b ~ -3.5°</td><td>12.4 ± 3.2</td><td>11.7 ± 2.6</td><td>17.7 ± 1.2</td><td>19.1 ± 1.5</td></tr><tr><td>b ~ +4.5°</td><td>8.6 ± 1.3</td><td>9.6 ± 1.1</td><td>8.3 ± 1.0</td><td>9.2 ± 1.1</td></tr><tr><td>b ~ -6°</td><td>7.8 ± 1.0</td><td>8.4 ± 0.8</td><td>10.6 ± 1.0</td><td>10.8 ± 1.1</td></tr><tr><td>b ~ -8.5°</td><td>7.0 ± 0.7</td><td>7.4 ± 0.6</td><td>8.7 ± 0.8</td><td>8.3 ± 0.8</td></tr></tbody></table>

is imprinted in the present data as well. Figure 13 shows the galactocentric radial velocity for all the stars, as a function of longitude, at fixed latitudes. Metal poor stars ([Fe/H] < -0.1) are shown in blue, metal rich stars ([Fe/H] > +0.1) in red. A linear relation, including the errors in radial velocity, was fit to the data in each panel, for the two metallicity components individually. Metal poor stars show a marginally slower rotation, that is, a flatter slope. Adopting the metallicity cuts mentioned above, and used elsewhere in this paper, the slopes are always compatible within the errors. However, if we move the metal poor cut at -0.3 dex, and the metal rich one at +0.2, that is, we exclude a larger metallicity bin as cross-contaminated, then the difference between the two slopes becomes more significant (Table 2). In other words, although the signal is marginally significant in these

<figure><img src="image_12.png" /><figcaption>**Fig. 11.** Velocity dispersion as a function of longitude, for fixed latitudes shown in different colors. Metal poor stars are shown in the left panel and metal rich ones in the right panel.</figcaption></figure>

data, the GIBS sample confirms the findings by Kunder et al. (2016) that metal rich stars rotate faster than metal poor ones.

## 6. Conclusions

We present the metallicity of approximately 5500 RC stars in 26 fields across the Galactic bulge. In all but two fields, the metallicity distribution is clearly bimodal with a metal poor component centered at [Fe/H]$_{MP}$ ≈ -0.4 dex and a metal rich component centered at [Fe/H]$_{MR}$ ≈ +0.3 dex. A small (few %)
<table><caption>**Table 3.** Components of the MDF fits and Nr of RC stars in each field.</caption><thead><tr><th rowspan="2">Field name</th><th colspan="3">Halo</th><th colspan="3">Metal poor</th><th colspan="3">Metal rich</th></tr><tr><th>σ</th><th>N/N<sub>tot</sub></th><th>σ</th><th>N/N<sub>tot</sub></th><th>σ</th><th>N/N<sub>tot</sub></th><th>σ</th><th>N/N<sub>tot</sub></th></tr></thead><tbody><tr><td>LRp8p4</td><td>-1.38</td><td>0.08</td><td>0.02</td><td>-0.30</td><td>0.31</td><td>0.54</td><td>0.24</td><td>0.17</td><td>0.44</td></tr><tr><td>LRp5p4</td><td>-</td><td>-</td><td>0.00</td><td>0.04</td><td>0.34</td><td>0.27</td><td>0.49</td><td>0.09</td><td>0.73</td></tr><tr><td>LRp0p4</td><td>-1.14</td><td>0.08</td><td>0.02</td><td>-0.31</td><td>0.27</td><td>0.56</td><td>0.19</td><td>0.18</td><td>0.42</td></tr><tr><td>LRm3p4</td><td>-1.05</td><td>0.15</td><td>0.03</td><td>-0.35</td><td>0.23</td><td>0.34</td><td>0.24</td><td>0.24</td><td>0.63</td></tr><tr><td>LRm8p4</td><td>-0.91</td><td>0.13</td><td>0.02</td><td>-0.47</td><td>0.20</td><td>0.47</td><td>0.16</td><td>0.21</td><td>0.50</td></tr><tr><td>LRp0m1</td><td>-1.83</td><td>0.25</td><td>0.00</td><td>-0.35</td><td>0.28</td><td>0.53</td><td>0.22</td><td>0.18</td><td>0.46</td></tr><tr><td>LRp8m2</td><td>-</td><td>-</td><td>0.00</td><td>-0.87</td><td>0.12</td><td>0.03</td><td>0.00</td><td>0.34</td><td>0.97</td></tr><tr><td>LRp3m2</td><td>-1.40</td><td>0.16</td><td>0.01</td><td>-0.23</td><td>0.24</td><td>0.43</td><td>0.36</td><td>0.18</td><td>0.56</td></tr><tr><td>LRp0m2</td><td>-0.95</td><td>0.13</td><td>0.04</td><td>-0.26</td><td>0.29</td><td>0.56</td><td>0.37</td><td>0.19</td><td>0.40</td></tr><tr><td>LRm5m2</td><td>-</td><td>-</td><td>0.00</td><td>-0.30</td><td>0.20</td><td>0.35</td><td>0.24</td><td>0.17</td><td>0.65</td></tr><tr><td>LRm8m2</td><td>-1.42</td><td>0.12</td><td>0.01</td><td>-0.29</td><td>0.24</td><td>0.43</td><td>0.17</td><td>0.17</td><td>0.56</td></tr><tr><td>HRp8m3</td><td>-</td><td>-</td><td>0.00</td><td>-0.45</td><td>0.21</td><td>0.17</td><td>0.24</td><td>0.19</td><td>0.83</td></tr><tr><td>HRp4m3</td><td>-0.89</td><td>0.26</td><td>0.05</td><td>-0.18</td><td>0.19</td><td>0.41</td><td>0.32</td><td>0.15</td><td>0.54</td></tr><tr><td>BW</td><td>-</td><td>-</td><td>0.00</td><td>-0.23</td><td>0.39</td><td>0.48</td><td>0.33</td><td>0.16</td><td>0.52</td></tr><tr><td>HRm5m3</td><td>-1.60</td><td>0.06</td><td>0.01</td><td>-0.39</td><td>0.25</td><td>0.16</td><td>0.21</td><td>0.16</td><td>0.83</td></tr><tr><td>HRm7m4</td><td>-1.53</td><td>0.05</td><td>0.01</td><td>-0.45</td><td>0.16</td><td>0.29</td><td>0.22</td><td>0.21</td><td>0.70</td></tr><tr><td>LRp8m6</td><td>-1.27</td><td>0.30</td><td>0.03</td><td>-0.27</td><td>0.21</td><td>0.52</td><td>0.29</td><td>0.20</td><td>0.45</td></tr><tr><td>LRp4m6</td><td>-0.82</td><td>0.09</td><td>0.03</td><td>-0.26</td><td>0.22</td><td>0.40</td><td>0.28</td><td>0.22</td><td>0.57</td></tr><tr><td>bulge6</td><td>-0.90</td><td>0.31</td><td>0.09</td><td>-0.35</td><td>0.17</td><td>0.41</td><td>0.08</td><td>0.16</td><td>0.50</td></tr><tr><td>LRm4m6</td><td>-0.77</td><td>0.15</td><td>0.06</td><td>-0.35</td><td>0.17</td><td>0.41</td><td>0.23</td><td>0.22</td><td>0.52</td></tr><tr><td>LRm8m6</td><td>-</td><td>-</td><td>0.00</td><td>-0.36</td><td>0.36</td><td>0.94</td><td>0.26</td><td>0.17</td><td>0.06</td></tr><tr><td>LRp8m8</td><td>-1.00</td><td>0.06</td><td>0.02</td><td>-0.34</td><td>0.24</td><td>0.80</td><td>0.32</td><td>0.09</td><td>0.18</td></tr><tr><td>LRp4m8</td><td>-1.25</td><td>0.09</td><td>0.01</td><td>-0.35</td><td>0.26</td><td>0.59</td><td>0.31</td><td>0.17</td><td>0.41</td></tr><tr><td>LRp0m8</td><td>-0.85</td><td>0.54</td><td>0.04</td><td>-0.37</td><td>0.22</td><td>0.81</td><td>0.33</td><td>0.17</td><td>0.16</td></tr><tr><td>LRm3m8</td><td>-1.01</td><td>0.13</td><td>0.05</td><td>-0.35</td><td>0.24</td><td>0.82</td><td>0.36</td><td>0.11</td><td>0.13</td></tr><tr><td>LRm8m8</td><td>-1.01</td><td>0.18</td><td>0.07</td><td>-0.41</td><td>0.18</td><td>0.65</td><td>0.13</td><td>0.27</td><td>0.28</td></tr><tr><td colspan="10">Global MDFs</td></tr><tr><td>b = -1°</td><td>-1.83</td><td>0.25</td><td>&lt;0.01</td><td>-0.35</td><td>0.28</td><td>0.53</td><td>0.22</td><td>0.18</td><td>0.46</td></tr><tr><td>b = -2°</td><td>-0.80</td><td>0.16</td><td>&lt;0.02</td><td>-0.24</td><td>0.27</td><td>0.49</td><td>0.28</td><td>0.21</td><td>0.49</td></tr><tr><td>b = -3.5°</td><td>-1.44</td><td>0.15</td><td>&lt;0.01</td><td>-0.31</td><td>0.31</td><td>0.33</td><td>0.26</td><td>0.20</td><td>0.67</td></tr><tr><td>b = +4°</td><td>-1.12</td><td>0.20</td><td>&lt;0.02</td><td>-0.35</td><td>0.27</td><td>0.45</td><td>0.24</td><td>0.22</td><td>0.54</td></tr><tr><td>b = -6°</td><td>-0.87</td><td>0.36</td><td>&lt;0.06</td><td>-0.32</td><td>0.25</td><td>0.51</td><td>0.22</td><td>0.23</td><td>0.43</td></tr><tr><td>b = -8.5°</td><td>-0.95</td><td>0.16</td><td>&lt;0.04</td><td>-0.38</td><td>0.22</td><td>0.74</td><td>0.30</td><td>0.19</td><td>0.23</td></tr></tbody></table>

population of stars with $[Fe/H] < -0.8$ is also identified in most of the fields, tentatively interpreted as halo contamination. The latter population is highly incomplete because our target selection box does not extend to the blue of the RC, therefore we did not attempt to characterize it.

Figure 14 shows the metallicity versus galactocentric radial velocity for all the stars in the present sample. Two populations in metallicity are clearly seen here. No stars were measured outside this figure's box, that is, we do not detect high velocity stars, inconsistent with the tail of the main velocity distribution. The two main bulge components peak at approximately the same metallicity in different fields, but the relative fraction of stars in each of them changes dramatically across the bulge area, giving rise to well-known radial gradients in the bulge mean metallicity. In the outer region studied here ($b = -8.5°$), the metal poor component largely outnumbers the metal rich one, with the former accounting for 73% of the total number of targets observed at this latitude, and the latter for the remaining 27%. The relative fraction of metal poor stars drops to 33% of the total

at $b = -3.5°$: a result in very good agreement with previous findings by the ARGOS (Ness et al. 2013a) and the Gaia-ESO (Rojas-Arriagada et al. 2014) surveys.

The present dataset allows us to reach latitudes closer to the Galactic plane than previous studies probing the population of the inner bulge. A first, important result of the present work is that the relative fraction of metal poor stars increases again from $b = -3.5°$ to $b = -2°$ and $b = -1°$, reaching 49% and 53% of the total in the latter two fields, respectively.

Having demonstrated that the bulge includes two metallicity components whose relative fraction changes across the bulge area, that is, they have a different spatial distribution, we mapped their density distribution by coupling the relative fractions estimated here with the bulge stellar density map derived in Valenti et al. (2016). The resulting density maps, shown in Fig. 9, demonstrate that the metal poor component has a spheroid-like spatial distribution, versus a “boxy” distribution of the metal rich component, and that the radial density gradient of the metal poor one is steeper. Due to the coarse spatial sampling
<figure><img src="image_13.png" /></figure>

**Fig. 12.** The trend of velocity dispersion with longitude, in six strips of constant latitude. The horizontal dashed line at $\sigma_{RV} = 100$ has been drawn to guide the eye. Notice that, in the outer bulge, $|b| > 3.5$, the metal poor component has a higher velocity dispersion compared to the metal rich one. However, at $b = -2$, $\sigma_{RV}$ is roughly equal for both metallicity bins, while at $b = -1$ the $\sigma_{RV}$ for the metal poor sample lies below that of the metal rich population. Notice that the metal poor component at $b = -3.5$ contains very few stars, for this reason the trend is noisy (dashed line).

of the GIBS fields, these maps have been derived imposing four-fold symmetry (about the Galactic plane, and about the bulge projected minor axis). Therefore, we are by construction unable to identify possible asymmetries related to the fact that at positive/negative longitudes we see the near/far side of the bar, respectively.

The two bulge metallicity components not only have different mean metallicity and different spatial distribution, they also have different kinematics. As already found by the BRAVA (Kunder et al. 2012) and the ARGOS survey (Ness et al. 2013b), in the outer bulge ($|b| > 4°$) the metal poor component has a higher radial velocity dispersion compared to the metal rich one, at all longitudes. A new result of the present investigation is that this behavior is reversed in the inner bulge. Specifically, at $b = -3.5°$ and $b = -2°$, the velocity dispersion of the metal poor stars becomes very similar to that of the metal rich ones, and at $b = -1°$, it becomes smaller. This is due to the fact that the velocity distribution of the metal poor component has a mild gradient with latitude ranging from $\sigma \sim 80 \text{ km s}^{-1}$ at $b = -8.5°$ to $\sigma \sim 125 \text{ km s}^{-1}$ at $b = -1°$. On the contrary, the metal rich component has a low velocity dispersion in the outer bulge ($\sigma \sim 60 \text{ km s}^{-1}$ at $b = -8.5°$) that increases to $\sigma \sim 145 \text{ km s}^{-1}$ at $b = -1°$. In other words, the central *sigma-peak* discovered in Paper I, is mostly due to the metal rich stars. The origin of this *sigma-peak* in the inner ~250 pc of the Galactic bulge has been investigated by Valenti et al. (2016), who demonstrated that it coincides with a peak in the stellar density (traced by counting the number of RC stars across the VVV survey area) and thus with a peak in the mass density.

Finally, we investigated whether the GIBS data show any evidence for the metal poor component rotating slower than

<figure><img src="image_14.png" /></figure>

**Fig. 13.** Galactocentric radial velocity versus longitude, at approximately fixed latitudes, as indicated in the labels. Metal poor stars ($[Fe/H] < -0.1$) are shown in blue and metal rich stars ($[Fe/H] > +0.1$) in red. A linear relationship has been fit in each panel separately for each metallicity component. Metal poor stars rotate marginally slower than metal rich ones.

<figure><img src="image_15.png" /></figure>

**Fig. 14.** Metallicity versus galactocentric radial velocity for the whole sample discussed here. The two metallicity components, that we refer to as metal rich and metal poor, respectively, are clearly seen here.

the metal rich one, as suggested by Kunder et al. (2016), who traced the former by means of RR Lyrae stars. Indeed, the metal poor component shows a marginally slower rotation also in the present dataset, although the difference is not compelling here.

Summarizing, we unambiguously detect two components in the Galactic bulge, each having a relatively narrow metallicity distribution, so that the global bulge MDF is clearly bimodal. These two components have a different spatial distribution, with the metal poor being more centrally concentrated, with a rounder density distribution as projected in the sky. On the contrary, the metal rich component is clearly “boxy”, as expected for a bar seen edge on. We stress that, although the present data are consistent with a spheroidal distribution of the metal poor component (although we did not investigate the 3rd component along the line of sight), we are not presenting evidence in favor of
it having a different origin (i.e., a classical bulge). Different models of bulge formation are being developed, arguing that a spheroidal shape can be obtained through different formation scenarios (Debattista, priv. comm.). Here we provide new observational constraints not previously available, that need to be fulfilled by a Milky Way bulge formation model.

**Acknowledgements.** We gratefully acknowledge support by the Ministry of Economy, Development, and Tourism's Millennium Science Initiative through grant IC120009, awarded to The Millennium Institute of Astrophysics (MAS), by Fondecyt Regular 1150345 and by the BASAL-CATA Center for Astrophysics and Associated Technologies PFB-06. This work was finalized at the Aspen Center for Physics, which is supported by National Science Foundation grant PHY-1066293. M.Z. and D.M. were partially supported by a grant from the Simons Foundation, during their stay in Aspen.

## References

Alard, C. 2001, A&A, 379, L44

Alonso, A., Arribas, S., & Martínez-Roger, C. 1999, A&AS, 140, 261

Alves-Brito, A., Meléndez, J., Asplund, M., Ramírez, I., & Yong, D. 2010, A&A, 513, A35

Athanassoula, E. 2005, MNRAS, 358, 1477

Babusiaux, C. 2016, PASA, 33, e026

Babusiaux, C., Katz, D., Hill, V., et al. 2014, A&A, 563, A15

Dékány, I., Minniti, D., Catelan, M., et al. 2013, ApJ, 776, L19

Gerhard, O., & Martinez-Valpuesta, I. 2012, ApJ, 744, L8

Gonzalez, O. A., Rejkuba, M., Minniti, D., et al. 2011a, A&A, 534, L14

Gonzalez, O. A., Rejkuba, M., Zoccali, M., et al. 2011b, A&A, 530, A54

Gonzalez, O. A., Rejkuba, M., Zoccali, M., Valenti, E., & Minniti, D. 2011c, A&A, 534, A3

Gonzalez, O. A., Rejkuba, M., Zoccali, M., et al. 2012, A&A, 543, A13

Gonzalez, O. A., Rejkuba, M., Zoccali, M., et al. 2013, A&A, 552, A110

Gonzalez, O. A., Zoccali, M., Vásquez, S., et al. 2015, A&A, 584, A46 (Paper II)

Hartwick, F. D. A. 1976, ApJ, 209, 418

Hill, V., Lecureur, A., Gomez, A., et al. 2011, A&A, 534, A80

Howard, C. D., Rich, R. M., Clarkson, W., et al. 2009, ApJ, 702, L153

Kunder, A., Koch, A., Rich, R. M., et al. 2012, AJ, 143, 57

Kunder, A., Rich, R. M., Koch, A., et al. 2016, ApJ, 821, L25

McWilliam, A. 2016, PASA, 33, e040

McWilliam, A., & Rich, R. M. 1994, ApJS, 91, 749

McWilliam, A., & Zoccali, M. 2010, ApJ, 724, 1491

Minniti, D., Olszewski, E. W., Liebert, J., et al. 1995, MNRAS, 277, 1293

Minniti, D., Lucas, P. W., Emerson, J. P., et al. 2010, New Astron., 15, 433

Mucciarelli, A., Pancino, E., Lovisi, L., Ferraro, F. R., & Lapenna, E. 2013, ApJ, 766, 78

Muratov, A. L., & Gnedin, O. Y. 2010, ApJ, 718, 1266

Nataf, D. M., Udalski, A., Gould, A., Fouqué, P., & Stanek, K. Z. 2010, ApJ, 721, L28

Nataf, D. M., Udalski, A., Skowron, J., et al. 2015, MNRAS, 447, 1535

Ness, M., & Freeman, K. 2016, PASA, 33, e022

Ness, M., Freeman, K., Athanassoula, E., et al. 2012, ApJ, 756, 22

Ness, M., Freeman, K., Athanassoula, E., et al. 2013a, MNRAS, 430, 836

Ness, M., Freeman, K., Athanassoula, E., et al. 2013b, MNRAS, 432, 2092

Ortolani, S., Renzini, A., Gilmozzi, R., et al. 1995, Nature, 377, 701

Pagel, B. E. J. 1989, in Evolutionary Phenomena in Galaxies, eds. J. E. Beckman, & B. E. J. Pagel, 201

Pasquini, L., Avila, G., Blecha, A., et al. 2002, The Messenger, 110, 1

Patsis, P. A., Skokos, C., & Athanassoula, E. 2002, MNRAS, 337, 578

Pietrukowicz, P., Kozłowski, S., Skowron, J., et al. 2015, ApJ, 811, 113

Ramírez, I., & Meléndez, J. 2005, ApJ, 626, 465

Rich, R. M. 1988, AJ, 95, 828

Rich, R. M., Origlia, L., & Valenti, E. 2007, ApJ, 665, L19

Rich, R. M., Origlia, L., & Valenti, E. 2012, ApJ, 746, 59

Rojas-Arriagada, A., Recio-Blanco, A., Hill, V., et al. 2014, A&A, 569, A103

Saha, K., & Gerhard, O. 2013, MNRAS, 430, 2039

Saito, R. K., Zoccali, M., McWilliam, A., et al. 2011, AJ, 142, 76

Schultheis, M., Cunha, K., Zasowski, G., et al. 2015, A&A, 584, A45

Stanek, K. Z., Mateo, M., Udalski, A., et al. 1994, ApJ, 429, L73

Terndrup, D. M. 1988, AJ, 96, 884

Udalski, A., Szymanski, M., Kubiak, M., et al. 2002, Acta Astron., 52, 217

Valenti, E., Zoccali, M., Gonzalez, O. A., et al. 2016, A&A, 587, L6

Vásquez, S., Zoccali, M., Hill, V., et al. 2013, A&A, 555, A91

Vásquez, S., Zoccali, M., Hill, V., et al. 2015, A&A, 580, A121

Wegg, C., & Gerhard, O. 2013, MNRAS, 435, 1874

Wegg, C., Gerhard, O., & Portail, M. 2015, MNRAS, 450, 4050

Zoccali, M., & Valenti, E. 2016, PASA, 33, e025

Zoccali, M., Renzini, A., Ortolani, S., et al. 2003, A&A, 399, 931

Zoccali, M., Hill, V., Lecureur, A., et al. 2008, A&A, 486, 177

Zoccali, M., Gonzalez, O. A., Vasquez, S., et al. 2014, A&A, 562, A66 (Paper I)