© ESO, 2012

1. Introduction

Since the seminal work of Bahcall & Soneira (1980a,b, 1984), the method of simulating star counts based on analytic Milky Way models has evolved considerably. The structural and evolutionary parameters of the Milky Way components have been refined in multi-population scenarios. However, a “concordance” model including the detailed shape and parameters of the density profiles, star formation history (SFH), age-velocity dispersion relation (AVR), etc., has not yet been reached.

On the other hand, there is a rapidly increasing set of observations consisting of more accurate and much deeper photometry to constrain the parameters of the models of each population. The observations include wider and deeper sky coverage and additional filter systems providing results that may require modifications of the existing detailed Milky Way models. The Sloan Digital Sky Survey (SDSS, York et al. 2000; Stoughton et al. 2002) covers more than one third of the celestial sphere providing photometric data in five bands of 2.6 × 10⁸ stars and spectra of more than 500 000 stars (for the most recent eighth data release (DR8), see Aihara et al. 2011). The SDSS has become one of the most widely used databases for Galactic and extragalactic astronomy in recent years. The availability of SDSS data enables us to make significant progress in constraining Galactic model parameters.

Analytic models are a powerful tool to constrain either evolutionary scenarios or assumptions about Galactic structure by comparing model predictions with a large variety of observational constraints such as star counts and kinematics. Models can be used to create probability distributions in the parameter space of observable or mock catalogues via Monte Carlo simulations. Realistic parameters and details for further improvement can be obtained by comparing model predictions with suitable observations, such as multi-directional photometric star counts. Within these parameters and functions, the density profiles, SFHs, and AVRs are important inputs.

To simulate and predict star counts in every direction of the sky in the bands of various photometric systems, several automated tools are available on the web, e.g., the “Besançon” model¹ by Robin et al. (2003) and the TRI-dimensional modeL of thE GALaxy (“TRILEGAL”) model² (Girardi & Bertelli 1998; Girardi et al. 2002, 2005), which are able to deal with isochrone construction, synthetic photometry, and simulations of the stellar populations of the main Galactic components. Sharma et al. (2011) presented Galaxia, a new tool which combines the advantages of the Besançon and the TRILEGAL model. The code generates very efficiently synthetic surveys and allows one to choose between a much larger variety of input models. The first version of Galaxia is expected to be presented online soon³.

The “Besançon model” (Robin et al. 2003) is presently the most elaborate and refined publicly available tool for predicting star counts. In this model, the Milky Way is divided into four components (thin disc, thick disc, spheroid, bulge), which are described by their spatial distribution, SFHs, initial mass functions (IMFs), sets of evolutionary tracks, kinematics, and metallicity characteristics, and include giants and a white dwarf population. From Monte Carlo simulations, mock catalogues including observable parameters as well as theoretical ones are obtained.

The TRILEGAL model is a Monte Carlo simulation to predict the probability distributions of stars across small sky areas and for special passbands systems, is still work in progress. It is a population synthesis code to simulate photometry of any field in the Milky Way (Girardi et al. 2005). It was developed to model population synthesis star counts of the Milky Way. It allows its users to create a pseudo stellar catalogue including positions, multi-colour photometry, and other physical parameters according to user-defined structural parameters.

On the basis of SDSS photometry from Data Release 6 (DR6, Adelman-McCarthy et al. 2008), Jurić et al. (2008) derived a three-dimensional (3D) stellar distribution of the Milky Way by applying universal photometric parallaxes. They then fitted two discs, each with exponential profiles in the radial and vertical directions. They found vertical scale heights of 300 pc and 900 pc for the thin and thick disc, respectively. The local density normalisation of the discs and basic parameters of the stellar halo were also determined. Subtracting the resulting, smooth Milky Way model, overdensities (and underdensities) were revealed that could be investigated in more detail.

We have developed a new local Milky Way model based on local stellar kinematics and star counts in the solar neighbourhood (Just & Jahreiß 2010, hereafter “Just-Jahreiß model”). Compared to the Besançon and TRILEGAL models, the main advantages of our approach are that vertical density profiles are fully consistent with the SFH and AVR in the gravitational field and that we use a well-calibrated empirical main-sequence (MS). In the current state of the Just-Jahreiß model, only MS stars are used instead of including the full Hertzsprung-Russell diagram for each stellar sub-population. The comparison of the basic Just-Jahreiß disc model with SDSS star counts in the NGP field with b > 80° enabled us to select a best-fitting SFH and AVR (Just et al. 2011). Predicted star-number densities fit the data with a median relative deviation of 5.6% across the Hess diagram. We use the Just-Jahreiß model to judge the ability of the TRILEGAL and the Besançon model to reproduce SDSS data at the NGP.

In this article, we present a test of the ability of the TRILEGAL code to reproduce the different Milky Way models mentioned above, since the TRILEGAL model is presently the only code that allows an interactive analysis of the impact of model parameter variations on star count predictions. Additionally, we try to determine a best-fit model by a systematic input-parameter optimisation in the TRILEGAL model. The Besançon model has no free parameters and we perform a multi-colour analysis over all five filters in the ugriz system. We restrict the analysis to the SDSS star-count data of the NGP field, since these are very high quality data and they give a direct measure of the controversial vertical structures of thin and thick discs.

The output catalogues of the TRILEGAL and Besançon models with photometric parameters can be used to analyse the contributions of different populations to the observed star counts and to quantify systematic deviations in order to identify the weaknesses that need to be resolved.

Section 2 introduces the data used in this work. In Sect. 3, we describe the details of the Just-Jahreiß, TRILEGAL, Besançon, and the Jurić models used in our analysis. In Sect. 4, we present the results of the TRILEGAL model concerning the isochrones, luminosity functions, and Hess diagrams. In Sect. 5, we present the comparisons between the Besançon model and the data using similar methods. In Sect. 6, we summarise the results of our study.

2. Data

The SDSS has become one of the most widely used databases for Galactic and extragalactic astronomy in recent years. The effective exposure time for imaging is approximately 54 seconds with a limiting magnitude of r ~ 22.6 mag. The g and r band point-source imaging of the SDSS ranges from 14 mag to 22.5 mag. We selected stars from the Catalog Archive Server (CAS) of the seventh data release of the SDSS (DR7, Abazajian et al. 2009). The SDSS table Star, which we utilise, contains the photometric parameters (no redshifts nor spectroscopic parameters) of all primary point-like objects classified as stars from PhotoPrimary. We only consider the NGP field (Galactic latitude b > 80°) in this publication. A series of lower latitude fields will be analysed in the future.

We use de-reddened magnitudes from the SDSS point-spread function photometry based on the values from the extinction map of Schlegel et al. (1998). In the following, all magnitudes and colours are de-reddened. We select all stars in DR7 within a colour range of (g − r) =  [−0.2,1.2]  and magnitude limits of g =  [14,20.5]  with magnitude errors below 0.2 mag in g and in r. In these magnitude and colour ranges, the error limits remove 0.88% of the sample in the NGP field. A total number of N = 276   180 stars remain in the final sample.

For our analysis, we produce Hess diagrams of the stellar number-density distribution in colour–magnitude bins normalised to 1 deg² sky area in the colour–magnitude diagram (CMD) for a given direction. One obvious deficiency are the unknown stellar distances, hence the vertical axis of the CMDs necessarily shows apparent magnitudes instead of absolute magnitudes. The combined distribution of (apparent) stellar luminosity and colour shows the spatial structure and the population properties in the chosen direction, which can then be compared with model-generated Hess diagrams.

In the top panel of Fig. 1, the Hess diagram in (g − r,g) is shown. Number densities per square degree, 1 mag in g and 0.1 mag in (g − r), are colour–coded on a logarithmic scale. The bin size is (0.05,0.5) in (g − r,g) with additional smoothing in steps of 0.01 mag. From the upper right (purple/blue regions) to the lower left (red triangular region) the thin disc, thick disc, and halo dominate, respectively (see Just et al. 2011). In the top left of the figure, no stars are visible. In the bottom left, white dwarfs, blue horizonal branch stars, and blue stragglers of the thick disc and halo may appear at (g − r) < 0.2, but the majority of objects are mis-identified point-like extragalactic sources.

The luminosity function is obtained by projecting the Hess diagram along the (g − r)-axis. In the lower panel of Fig. 1, the SDSS data are shown for the complete range in the g band in order to illustrate the total numbers of stars and the cutoffs at bright and faint magnitudes. Poissonian error bars were adopted. The incompleteness at g > 22 mag is obvious and a significant confusion with mis-identified extragalactic sources is likely to occur here (York et al. 2000). To be on the safe side, we choose for our analysis g = 20.5 mag as the faint limit. One significant jump appears at g = 14 mag, which is caused by the incompleteness due to the large magnitude errors. The bright limit at 14 mag is set by the saturation of the stars in the SDSS observations.

Apart from the DR7 catalogue data, we select the fiducial isochrones of five globular clusters determined by An et al. (2008) to compare the photometric properties of the thick disc and halo populations in the TRILEGAL and Besançon models. Using the extinction parameters and distance moduli of these globular clusters (Harris 1996), the de-reddened absolute magnitudes and colour indices are determined. The parameters of the globular clusters are listed in Table 1. These clusters are selected to cover the metallicity range of the thick disc and halo adopted in the TRILEGAL and the Besançon models.

3. Models

We describe the basic properties of the four models compared with the SDSS data.

3.1. Just-Jahreiß model

Just & Jahreiß (2010) presented a local disc model where the thin disc is composed of a continuous series of isothermal sub-populations characterised by their SFH and AVR. A distinct isothermal thick disc component can be added. The vertical density profiles ρ_s,j(z) of all sub-populations j are determined self-consistently in the gravitational potential Φ(z) of the thin disc, thick disc, gas component, and dark matter halo. The density profiles of the MS stars are calculated by adding the profiles of the sub-populations, accounting for their respective MS lifetimes. In a similar way, the thicknesses and the velocity distribution functions of MS stars are derived as a function of lifetime.

For the thin disc of the Just-Jahreiß model, a factor of 1/g_eff = 1.53 with g_eff = 0.654 needs to be applied in order to convert the present-day stellar density to the density of ever formed stars. For the thick disc and halo, the conversion factor is 1.8.

In Just & Jahreiß (2010), we showed that the vertical profile of an isothermal thick disc can be well-fitted by a sech^α_t(z/α_th_t) function. The power-law index α_t and the scale height h_t are connected via the velocity dispersion.

The SFH and the AVR are the main input functions that are to be determined as a pair of smooth continuous functions. For each given SFH, the AVR is optimised to fit the velocity distribution functions of MS stars in V-band magnitude bins in the solar neighbourhood. A metal enrichment law (i.e., the age-metallicity relation (AMR)) is included that reproduces the local G-dwarf metallicity distribution. Additionally, a consistent IMF is determined, that reproduces the local luminosity function. In Just & Jahreiß (2010), four different solutions with similar minimum χ² but very different SFHs were discussed.

In Just et al. (2011), we combined the models of Just & Jahreiß (2010) with the empirically determined photometric properties of the local MS in the ugriz filter system. By adopting a best-fitting procedure to determine the local normalisations of thin disc, thick disc, and stellar halo as a function of colour, we compared the star count predictions of the different models for the NGP field. In this procedure, the thick disc parameters α_t and h_t were optimized and we added a simple stellar halo with a flattened power-law-density distribution.

To measure the quality of the model and fit the data, we here define a χ² value for the Hess diagram. Each Hess diagram contains m × n boxes of magnitudes and colours. The relative difference (data − model)/model of each Hess diagram is a matrix D_ij with the indices i for colour and j for magnitude.

The total χ² is used to measure of the modelling quality, which is our standard way of evaluating the goodness of a model. We apply the same definition of χ² as Just et al. (2011) to measure the goodness of the TRILEGAL and Besançon models. We also define the median deviation Δ_med = median(| D_ij|) of the relative differences D_ij for each Hess diagram to quantify the typical difference between model and data.

In Just et al. (2011) it turned out that the Just & Jahreiß (2010) “model A” is clearly superior to the other models. We therefore select “model A” as the fiducial Just-Jahreiß model. Its basic density profiles and parameters are listed in the last column of Table 3 for reference.

3.2. TRILEGAL

The main goal of the TRILEGAL model is the simulation of Milky Way star counts in several photometric systems using simple profiles. For this purpose, Girardi et al. (2002) developed a detailed method to convert theoretical stellar spectral libraries to photometrically measurable luminosities. These authors can produce these transformations routinely for any new photometric system, such as the five-colour filters of the SDSS (Girardi et al. 2004). In addition, TRILEGAL has been used to predict stellar photometry in different systems.

Being able to simulate both very deep and very shallow surveys is one of the most difficult goals of TRILEGAL. For shallower surveys such as Hipparcos data, MS, red giant branch, and red clump stars with partly very low number statistics need to be generated. On the other hand, for more sensitive surveys an extended lower-luminosity MS, which reaches down to M_V ~ 30 mag corresponding to stellar masses of M ~ 0.1   M_⊙ as well as old white dwarfs need to be included.

The measurements of stellar age (τ) and metallicity (Z) provide parameters such as bolometric magnitude M_bol, effective temperature T_eff, surface gravity g, core mass, and surface chemical composition.

The IMF is a crucial ingredient in Galactic modelling, because it determines the relative number of stars with different photometric properties and lifetimes. For comparison and convenience, the default shape of the IMF is a log-normal function in TRILEGAL, following Chabrier (2001).

The SFH describes the age distribution of the stars. For old components, the shape of the SFH is quite poorly constrained, because the observed star counts are quite insensitive to it. For the thin disc, the SFH has not yet been well-determined (e.g., Aumer & Binney 2009; Just & Jahreiß 2010). The interface of TRILEGAL allows only the choice between a constant star formation rate (SFR) up to an age of 11 Gyr and a two-step SFH with a 1.5 times enhancement in the SFR between the ages of 1 and 4 Gyr.

The AMR is also an important ingredient of the Galactic model, because it can help modify the stellar CMDs. For the thin disc, the abundances of Rocha-Pinto et al. (2000) are adopted in TRILEGAL.

There are five Galactic components in TRILEGAL: thin disc, thick disc, halo, bulge, and an extinction layer in the disc.

The mass density of the thin disc decreases exponentially with Galactocentric radius R projected onto the plane of the disc. The Galactic thin disc is radially exponential and locally not isothermal, because the scale height h_d is a function of stellar age and the thin disc is not isolated. In TRILEGAL, one can choose between an exponential and a sech² profile for the vertical densities of both the thin and the thick disc. The typical height of the disc sub-populations depends on stellar age. The stellar density in the solar neighbourhood ρ_d(⊙ ) is the mass density of all stars ever born there.

The density of the thick disc can be described as either a double exponential or vertically with a sech² form similar to that of the thin disc. However h_t is not a function of stellar age because the thick disc is predominantly an old population. The SFH is constant over an age range of 11–12 Gyr with a metallicity [Fe/H]  = −0.67 ± 0.1 dex and an α-enhancement of [α/Fe]  ~ 0.3.

The mass density ρ_h of the halo follows a de Vaucouleurs r^1/4 law (de Vaucouleurs 1948) and allows for a spheroidal flattening (Gilmore 1984). The density profile of the stellar halo can be described by the parameters effective radius, oblateness, and local density (Young 1976). The SFH is constant over an age range of 12–13 Gyr and there are two choices for the metallicity distribution. We fix the metallicity to be [Fe/H]  = −1.6 ± 1 dex with a corresponding α-enhancement of [α/Fe] = 0.3.

The bulge has not yet been calibrated well in TRILEGAL because of a lack of observational constraints. Since the bulge does not contribute to high Galactic-latitude data, we do not discuss its parameters any further.

We do not use the dust extinction model of TRILEGAL but instead we compare star counts with de-reddened observational data as described in Sect. 2.

Each simulated photometric catalogue in TRILEGAL has a given direction that is a conical area with a maximum sky coverage of 10 deg². The centre of the area can be specified in either Galactic coordinates (l,b) or equatorial coordinates (α,δ). The resolution in magnitudes,  △ m, in TRILEGAL’s calculations can be modified by the user. The default value is  △ m = 0.1 mag and can be reduced to a minimum of 0.05 mag. Any data point with a resolution less than  △ m will not be presented. The output catalogues of TRILEGAL are based on a random number generator, which allows one to investigate the effect of noise.

Since the maximum sky area is much smaller than the observed field (313.36 deg² for the NGP field), we generate catalogues for 37 different sight lines (Table 2) to cover the NGP field. We add up all the star counts and then normalise to star counts per deg².

3.3. Besançon model

Robin et al. (2003) developed the population synthesis approach to simulate the structure and evolution of the Milky Way. The Besançon model is based on different well-calibrated data sets and the structural parameters are optimised to reproduce a set of selected fields on the sky. The construction of each component is based on a set of evolutionary tracks and assumptions about stellar number-density distributions, which are both constrained by either dynamical considerations or empirical data and guided by a scenario of a continuous formation and evolution of stellar populations.

The Besançon model deals with the thin disc is composed of seven isothermal sub-populations with different age ranges, kinematics parameters (W velocity dispersion σ_W), scale heights, and local densities. The thick disc is considered as a single population with an age of 11 Gyr and a local density value (1.34 × 10^-3   M_⊙   pc^-3). The halo is modelled by a single population of older age (14 Gyr) and lower local density (2.2 × 10^-5   M_⊙   pc^-3 including white dwarfs). This prescription is similar to our Just-Jahreiß model (Just & Jahreiß 2010). The main differences are that 1) in the Just-Jahreiß model, we adopt a continuous series of isothermal sub-populations in the thin disc; 2) the density profile and local density of each component are determined by solving iteratively Poisson and Jeans equations that are consistently based on the local kinematic data and the SDSS photometric data; and 3) we use a mean empirically calibrated MS instead of a full Hertzsprung-Russell diagram (HRD) based on stellar evolution.

Since we limit our targets to the NGP field with a Galactic latitude b > 80°, stars belonging to the Galactic bulge can be ignored. The stellar populations in the NGP field can be divided into three distinct components: the thin disc, the thick disc, and stellar halo. For details of the properties of the Galactic components, we refer to Robin et al. (2003). With the NGP field, we can test the vertical structure of the Galactic components and the calibration of the underlying isochrones in the SDSS filter system.

The Besançon model outputs u^∗g′r′i′z′ photometry. To obtain SDSS magnitudes in g, r, i, and z and the colours u − g, g − r, r − i and i − z from the Besançon model, we transform the photometric system from the u^∗, g′, r′, i′ and z′ filter set using a two-step method. In the first step, we transform from the special u^∗ filter to u′ (Clem et al. 2008). In the second step, we transform from the u′, g′, r′, i′, and z′ filters to u, g, r, i, and z following the transformations in Tucker et al. (2006) and Just & Jahreiß (2008).

3.4. Jurić model

In Jurić et al. (2008), a map of the 3D number density distribution in the Milky Way derived by photometric parallaxes that were based on star count data of SDSS DR6 is presented. The stellar distances in the sample were determined by the colour–absolute magnitude relation of MS stars. Most stars are at high latitudes and cover a distance range from 100 pc to 20 kpc and 6500 deg² on the sky in the sample.

The density distribution is fitted using two double-exponential discs and the sample implies the existence of an oblate halo. The estimated errors in the model parameters are smaller than  ~20% and the errors in the scale heights of discs are smaller than  ~10%. By adjusting the fraction of binaries and the bias correction, the scale heights of the discs are found to be h_d = 300 pc and h_t = 900 pc, and the local calibration is ρ_t(⊙ )/ρ_d(⊙ ) = 0.12. The best-fit of the spheroid leads to an oblateness of 0.64, a power-law radial profile of ρ_h ∝ r^-2.8, and a local calibration of ρ_h(⊙ )/ρ_d(⊙ ) = 0.005. After subtraction of the smooth background distribution, there are still overdensity areas, in particular areas that are explained as tidal streams including the “Virgo overdensity” (Duffau et al. 2006) and the Monoceros stream (Yanny et al. 2003).

The solar offset Z_⊙ is (25 ± 5) pc or (24 ± 5) pc at the bright and faint ends, respectively, which is consistent with the value of 24.2 pc adopted by TRILEGAL.

4. Simulations and results of TRILEGAL

The model TRILEGAL can be accessed via a web-based interface. Users can simulate the stellar populations in a special volume defined by the luminosity range and direction on the sky using different parameter sets. Because the parameter sets and profile shapes are relatively simple and the web-based programme can return the results in real time, we can compare the output interactively with the data.

In general, it is difficult to find the “best” input parameters for the components by using the web interface interactively, since the transitions of thin to thick disc and thick disc to halo in the CMD are not known a priori. We cannot search independently for the best-fitting parameters in each component. We therefore first test both the parameters of the default model and those that reproduce best the other three Milky Way models described in Sect. 3.

We fix some parameters and input functions in the TRILEGAL model in order to facilitate the comparison with observational data. We fix the IMF to a Chabrier log-normal function and use the default binary distribution with a binary fraction of 0.3 and a mass-ratio range of 0.7–1. The solar position is at the default value of (R_⊙,z_⊙) = (8700,24.2) pc and we fix the SFH for the thin disc to be constant. For the other input parameters in the TRILEGAL model, we selected five sets of parameters as given in Table 3. The “set 1” parameters in this table are from the default input of TRILEGAL in the web interface, “set 2” is best adjusted to the Just-Jahreiß model (Just & Jahreiß 2010), “set 3” to the Besançon model (Robin et al. 2003), “set 4” to the Jurić model (Jurić et al. 2008), and “set 5” is optimised (as described in Sect. 4.6).

Because in the TRILEGAL model, we have to choose between either an exponential or a sech² shape of the disc profiles $\begin{array}{ccc}{\mathit{\rho }}_{\mathrm{t}}\mathrm{\left(}\mathit{z}\mathrm{\right)}& \mathrm{=}& {\mathit{\rho }}_{\mathrm{t}}\mathrm{(}\mathrm{\odot }\mathrm{)}\exp \left(\mathrm{-}\frac{\mathit{z}}{{\mathit{h}}_{\mathrm{t}}}\right)\mathit{,}\\ {\mathit{\rho }}_{\mathrm{t}}\mathrm{\left(}\mathit{z}\mathrm{\right)}& \mathrm{=}& {\mathit{\rho }}_{\mathrm{t}}\mathrm{(}\mathrm{\odot }\mathrm{)}{\mathrm{sech}}^{\mathrm{2}}\left(\frac{\mathit{z}}{\mathrm{2}{\mathit{h}}_{\mathrm{t}}}\right)\mathit{,}\end{array}$we determine for all combinations the minimum χ² values and then select the profiles with the lowest χ² value.

At colours bluer than g − r = 0.2 mag, the SDSS data are unreliable owing to contamination by extragalactic objects and all selected models are unable to generate sufficient stars in this region. Therefore, the mean χ² value  ⟨ χ² ⟩  of the reduced colour range 0.2 < g − r < 1.2 is treated as the standard way of determining the quality of the model. For completeness, the values  ⟨ χ² ⟩  ^⋆  for the full colour range are also listed in Table 3.

Figure 2 compares the vertical density profiles of the thick disc and halo for the different sets. There are significant differences between the models, even in the distance regime where the thick disc component dominates the star counts in the Hess diagram. We discuss this point later.

To separate the impact of the synthetic colours and luminosities of the underlying stellar populations from the structural parameters of the Milky Way model, we compare the colour–absolute magnitude diagrams provided by TRILEGAL with the observed fiducial isochrones of globular clusters with corresponding metallicities. This comparison is only useful for both thick disc and halo, which are essentially simple stellar populations with small ranges in age and metallicity. The colour–absolute magnitude diagrams of Set 2 are plotted in Fig. 3. The cluster M 71 ([Fe/H]  = −0.73) is used to provide fiducial isochrones of the thick disc (see An et al. 2008). The three clusters M 71, M 13, and M 92 are used to provide fiducial isochrones of the halo whose mean metallicity is given as −1.6 (Girardi et al. 2005). The consistency between the model and the fiducial isochrones for the location of the MS and the turn-off implies that the stellar population synthesis of TRILEGAL is satisfactory for predicting the dominant contribution to the star counts. However, there are some systematic inconsistencies concerning the location of the giant branches of the stellar isochrones.

4.1. Default input (set 1)

The default input of TRILEGAL is adopted directly as parameter set 1. The parameters are listed in Col. 3 of Table 3. The default set 1 has no thick disc. For the vertical structure pointing to the NGP, set 1 simulates star counts using only thin disc and halo.

The top row of Fig. 4 shows the luminosity functions in the g and r bands of the SDSS data and of the models in the colour range 0.2 < g − r < 1.2. The default model (set 1) fits the data very well in the selected magnitude range 14–20.5 mag.

In the CMD, the picture is different. We construct a Hess diagram that is similar to the observed Hess diagram in Fig. 1. The top row of Fig. 5 shows the analysis of the Hess diagram. The left panel shows the simulated Hess diagram, which we compare directly to the observed Hess diagram. The middle panel shows the median absolute relative differences D_ij with a colour–coded range of  ± 0.5. The right panel gives the unreduced χ² distribution on a logarithmic scale. From the relative differences, some clear features can be identified (middle panel of the top row of Fig. 5):

• 1.

  The large white “island” in the middle of the diagram is a consequence of the lack of a thick disc in the model.

• 2.

  The missing stars at the top-right corner of g − r ~ 1.1 and g ~ 15 mag are a hint of an underestimation of the local thin-disc density caused by the strongly flattened sech² profile close to the Galactic plane.

• 3.

  The underestimated number counts at the faint blue end of 0.2 < g − r < 0.4 and g > 19.5 mag show that the halo density profile falls off too steeply.

• 4.

  The horizontally adjacent black-and-white coloured areas in the colour range 0.15 < g − r < 0.45 at g ~ 17 mag show that the (thick disc) population is replaced in the model by a more metal-poor (halo) population with a bluer F-turnoff region. This is additional strong evidence that a thick disc component of intermediate metallicity is missing in the model.

• 5.

  At the faint end (0.15 < g − r < 0.25 and 15 < g < 19), the stars missing in blue colours (halo) are compensated for by too many (or slightly too blue) M dwarfs of the thin disc in the luminosity function.

The value of the mean reduced  ⟨ χ² ⟩  = 151.98 is significantly larger than that of the Just-Jahreiß model (=4.31, Table 3). The median relative deviation of the data from the model in the Hess diagram Δ_med is 0.3035, which is 5.4 times larger than for the Just-Jahreiß model.

A comparison of Figs. 4 and 5 indicates that the analysis of the luminosity function integrated over a large colour range alone may result in misleading conclusions about the underlying stellar populations.

4.2. Just-Jahreiß input (set 2)

For parameter set 2, we derive the best representation of the input function of the Just-Jahreiß model with the TRILEGAL parameters. With this “model of a model”, we investigate the capability of TRILEGAL to reproduce a Milky Way model suggested by other authors.

4.2.1. Input parameters

For parameter set 2, which is based on the Just-Jahreiß model, we derived a best-fit TRILEGAL profile of the thick disc to the sech^αt profile. Both the exponential and the sech² profiles are fitted with a standard routine and the unreduced χ² values are used to select the best fit. The exponential form results in χ² = 4.0 × 10^-8, which is a better fit than the sech² profile with χ² = 6.6 × 10^-8. We therefore use the exponential profile. The scale height of the thick disc is 793 ± 6 pc and the local density is 0.0029 ± 0.000015 M_⊙   pc^-3. The derived local surface density of the thick disc is Σ_t(⊙ ) = 4.63 M_⊙   pc^-2.

In TRILEGAL, the scale height of the thin disc is a function of both age and the parameters z₀, t₀, and α (Girardi et al. 2005). The parameters are determined by fitting the numerical Just-Jahreiß model relationship of the thin disc. We fit the thin disc parameters in exponential and sech² form to the Just-Jahreiß model, respectively. The latter fit results in a smaller χ² (namely 4.7 × 10^-6) than the former one (8.3 × 10^-6), hence we choose the sech² form as the profile of the thin disc. The derived local surface density of the thin disc is Σ_d(⊙ ) = 32.26 M_⊙   pc^-2, which is 10% higher than the surface density of the thin disc of the Just-Jahreiß model.

The halo of the Just-Jahreiß model is fitted by the TRILEGAL form in the oblate r^1/4 spheroid. The best-fit parameters are local density ρ_h(⊙ ) = 0.001576 M_⊙   pc^-3, oblateness q_h = 0.65, and effective radius r_eff = 11994 pc.

To take into account the estimated total local density of all stars ever formed, the present-day densities of the Just-Jahreiß model were corrected accordingly (see Sect. 3.1).

4.2.2. Star counts

The second row in Figs. 4 and 5 shows the star count results of TRILEGAL based on the Just-Jahreiß parameters. The luminosity functions are of poorer quality than that for the default parameter set 1, because the halo is now over-represented by about 0.2 dex. In the CMD, the thick disc fills the minimum in the centre of the Hess diagram seen in the default model. The F-turnoff of the thick disc is still slightly too red (metal-rich). In the faint (distant) regime, there are significantly too many stars in all components.

The overall fit of the Hess diagram with  ⟨ χ² ⟩  = 214.65 is worse than for the default model owing to the poor fit in the faint magnitude range. The median absolute relative difference in Hess diagram Δ_med is 0.2993.

4.3. Besançon input (set 3)

4.3.1. Model parameters

The thick disc of Robin et al. (2003) is described by a parabolic shape up to |z| = 400 pc and an exponential profile beyond this height.

We fit this piecewise-defined function using a sech² and an exponential function, respectively. The exponential function with a smaller fitting error is adopted in Col. 5 of Table 3. The scale height of 910 pc and a local density of 0.0015 M_⊙   pc^-3yield a local surface density Σ(⊙ ) = 2.82 M_⊙   pc^-2.

There are seven sub-populations in the thin disc of the Besançon model. Each sub-population has its own age range, metallicity, and spatial vertical distribution. We fit scale heights and local densities of each sub-population (see Table 4). Exponential profiles yield smaller χ² values.

In the next step, the scale heights as a function of age are fitted by the TRILEGAL law. The best-fit parameters z₀, t₀, and α are 13.61 pc, 0.02016 Gyr, and 0.5524, respectively. The total local density ρ_d(⊙ ) is the sum of ρ_i(⊙ ) of each population, namely 0.030 M_⊙   pc^-3. The derived local surface density of the thin disc is 25.20 M_⊙   pc^-2. The density was corrected by the factor 1/g_eff = 1.53 to obtain the mass of ever formed stars.

The local density of the halo in the Besançon model is determined by the SFH and the IMF. The adopted value is much lower than our local calibration established using SDSS data. We therefore fit the halo of the Besançon model with the effective radius, oblateness, and local density of the TRILEGAL de Vaucouleur profile in two steps:

• 1.

  we fit the vertical profile to find the effective radius and the oblateness;

• 2.

  we use the local density as a new fitting parameter to reproduce the halo star counts of the SDSS data and then use this value to fit a model of the parameter set 3 (Table 3).

In the first step, we find the best-fit parameters q_h = 0.71 and r_eff = 39.436 kpc. In the second step, we determine the local density of the halo by fixing q_h and r_eff and fitting the SDSS data in a restricted area of the Hess diagram where the halo dominates (i.e., in the dashed triangle of Fig. 6, see Sect. 4.6). The best-fit local density of the stellar halo is ρ_h(⊙ ) = 0.00009 M_⊙   pc^-3, which is about ten times higher than the value used in the Besançon model.

4.3.2. Star counts

The third row in Figs. 4 and 5 shows the star count predictions of TRILEGAL based on the Besançon parameters. The luminosity functions appear to disagree significantly with these predictions at bright magnitudes primarily owing to the presences of the thin and thick discs, but be a reasonably fit at the faint end. In the CMD, the thin and thick disc regime are strongly under-represented in the model. The halo with the original local density also appears to be under-represented.

The overall fit to the Hess diagram with  ⟨ χ² ⟩  = 207.60 is tighter than that of parameter set 2 but worse than that of the default set (set 1). The Δ_med is 0.3445, which is worse than for sets 1 and 2.

This result is not a judgement on the Besançon model compared to the Just-Jahreiß or TRILEGAL model. It is a measure of the flexibility of the TRILEGAL model. A comparison of the Besançon model with the NGP data is later discussed in more detail.

4.4. Jurić input (set 4)

4.4.1. Model parameters

We do not consider the scale height as a function of the age of the subpopulation in the Jurić case (α = 0). Jurić et al. (2008) fit the basic parameters of their Galactic model via star count observations from the SDSS. Only the final parameters of their discs are considered here. Jurić et al. (2008) describe the vertical profiles of the thin and thick disc by simple exponentials with h_d = 300 pc and ρ_d(⊙ ) = 0.029 M_⊙   pc^-3 for the thin disc and h_t = 900 pc, ρ_t(⊙ ) = 0.0034 M_⊙   pc^-3 for the thick disc. This corresponds to 12% of the local density of the thin disc given by Jurić et al. (2008). The derived local surface density of the thick disc is Σ_t(⊙ ) = 6.12 M_⊙   pc^-2. The density has been corrected to the mass of ever formed stars with the factor g_eff for both the thin and thick discs.

The halo density profile is fitted by TRILEGAL with three parameters. The best-fit parameters are an oblateness of q_h = 0.59, a local density of ρ_h(⊙ ) = 0.000145, and an effective radius of r_eff = 18304 pc.

4.4.2. Star counts

The bottom rows in Figs. 4 and 5 show the star count results of TRILEGAL based on the Jurić parameters. The luminosity functions are slightly more clisely reproduced than for set 3. They show a significant deficiency at bright magnitudes dominated by the thin and thick disc and too many stars at the faint end. In the CMD, the fit of the thick disc regime and the halo is reasonable in the model. At the red end, the bright part is under-represented and the faint end is over-represented.

The overall fit of the Hess diagram with  ⟨ χ² ⟩  = 137.88 is the best of the four parameter sets but has the worst Δ_med (0.4639).

4.5. Discussion

An overall comparison of sets 1 to 4 shows that the default parameter set fits best the luminosity functions in the g and r band. On the one hand, this confirms the fitting procedure used in TRILEGAL to determine these parameters mainly by comparing luminosity functions in different fields. On the other hand, it is remarkable that the thick disc component, which becomes obvious in the Hess diagram analysis, is extremely poorly fit.

All other model adaptations fail to closely reproduce the luminosity functions. No model with a thin disc, thick disc, and halo can reproduce the shallow regime of the luminosity function in the range 17   mag < g < 19   mag. The analysis of star counts in the CMD shows that disagreements between the data and models exceed  ± 50% in some regions of the Hess diagram. In this sense, all models fail to provide a good fit of the CMD.

The default model shows that luminosity functions alone are inappropriate for determining the structure of the Milky Way owing to the large distance range in each component. Nevertheless, the total χ² value for the Hess diagram fit is relatively small compared to the other sets.

In the next section, we try to optimise the TRILEGAL parameters to find the best model for the CMD at the NGP.

4.6. Optimisation (set 5)

To find the best-fit TRILEGAL model, it is necessary to identify the best-fitting parameters. Because there are so many parameters that can be adjusted within TRILEGAL, we cannot determine them simultaneously. However since the thin disc, thick disc, and the halo do not completely overlap in the (g − r,g) CMD, it is possible to optimise the parameters for the halo, thick disc, and thin disc sequentially.

Figure 6 gives the schematic diagram for distinguishing areas where only the halo or both the halo and thick disc contribute significantly to the star counts, and where only the thick and thin discs contribute in this manner. We define the separating lines in terms of the relative densities of these three components based on the Just-Jahreiß model. In the resulting triangle (dashed line), the halo dominates over the other two components because the relative density of the discs to the halo is less than 1%. In the area below the solid line, the thin disc can be neglected, because the relative density of thin disc to thick disc is less than 1%. We fit first the halo parameters in the triangle. We then fix the halo parameters and fit the thick disc parameters below the full line by including the star counts of the halo. Finally, we investigate the thin disc parameters by taking into account all three components in the full CMD.

4.6.1. Halo fitting

Locating an area with halo stars but very few disc stars in the CMD is the first step. Our choice is marked by the dashed triangle in Fig. 6. The discs stars in this selected region can be neglected, i.e., (ρ_thin + ρ_thick)/ρ_halo < 0.01 as estimated by the Just-Jahreiß model. The boundary conditions to determine the halo parameters are g − r ≥ 0.2 and g ≥ 15.5 + 12.5(g − r). We choose an oblate r^1/4 spheroid and search for a best fit in the parameter ranges $\begin{array}{ccc}& & \mathrm{0.00005} {\mathit{M}}_{\mathrm{\odot }}\mathit{/}{\mathrm{pc}}^{\mathrm{3}}⩽{\mathrm{\Omega }}_{\mathrm{h}}\mathrm{\left(}\mathrm{\odot }\mathrm{\right)}⩽\mathrm{0.00035} {\mathit{M}}_{\mathrm{\odot }}\mathit{/}{\mathrm{pc}}^{\mathrm{3}}\mathit{,}\\ & & \mathrm{2000} \mathrm{pc}⩽{\mathit{r}}_{\mathrm{eff}}⩽\mathrm{30}\hspace{0.17em}\mathrm{000} \mathrm{pc}\mathit{,}\\ & & \mathrm{0.57}⩽{\mathit{q}}_{\mathrm{h}}⩽\mathrm{0.70.}\end{array}$\arraycolsep1.75ptThe fitting result is shown in Fig. 7. The best-fit halo parameters are listed in the last column (set 5) of Table 3 with  ⟨ χ² ⟩  = 4.86.

4.6.2. Thick disc fitting

We fix the parameters of the halo determined in Sect. 4.6.1 to fit the parameters of the thick disc in the area dominated by the thick disc and halo (below the solid line determined by g ≥ 16.5 + 3.64(g − r) in Fig. 6). In this area, only relatively few thin-disc stars exist, i.e., ρ_thin/(ρ_thick + ρ_halo) < 0.01.

The sech²() and exp() prescriptions are tried to fit the vertical profile of thick disc with the parameters scale height h_t and local calibration ρ_t(⊙ ) within the following ranges: $\begin{array}{ccc}& & \mathrm{0.0010} {\mathit{M}}_{\mathrm{\odot }}\mathit{/}{\mathrm{pc}}^{\mathrm{3}}⩽{\mathit{\rho }}_{\mathrm{t}}\mathrm{\left(}\mathrm{\odot }\mathrm{\right)}⩽\mathrm{0.00035} {\mathit{M}}_{\mathrm{\odot }}\mathit{/}{\mathrm{pc}}^{\mathrm{3}}\mathrm{;}\\ & & \mathrm{400} \mathrm{pc}⩽{\mathit{h}}_{\mathrm{t}}⩽\mathrm{900} \mathrm{pc}\mathit{.}\end{array}$\arraycolsep1.75ptThe minimum  ⟨ χ² ⟩  = 16.22 was found for the sech²() profile with the best-fitting parameters listed in Col. 7 of Table 3. The minimum  ⟨ χ² ⟩  of the fit using an exponential function is 0.312.

The fit result, shown in Fig. 8, shows on the one hand the homogeneous quality of the fit, but on the other hand a systematic trend along the MS of the thick disc with too few stars in the G-dwarf regime and too many K and M dwarfs.

4.6.3. Thin disc choice

For the thin disc fitting, we now fix the halo and thick disc as determined above. For the thin disc, we cannot get a best fit for all the free parameters. On the basis of the analysis of the advantages and disadvantages of the thin disc models of for the parameter sets 1, 2, 3, and 4, we adjust the parameters iteratively to get a better result. All iterations to improve the parameters in order to find an optimal Milky Way model must be done by hand. The power-law index α in the scale-height function changes the fraction of different sub-populations as a function of height above the mid-plane. In the Hess diagram, this affects the balance of densities between g − r > 1 (faint, nearby) and 0.2 < g − r < 0.4 (bright, more distant). The values for z₀ and ρ_d,0 affect the total density of the thin disc. These are the general directions for adjusting the parameters. The value of t₀ was not changed.

We performed these investigations for both vertical profile shapes. The sech² profile results in a better fit with a minimum  ⟨ χ² ⟩  = 1.33. (The exponential profile leads to a best fit with  ⟨ χ² ⟩  = 2.95).

4.6.4. Model properties

The Hess diagram analysis for input parameter set 5 is plotted in Fig. 9. The χ² distribution is smoother than for the other four sets, but the systematic features in the relative deviations of data from model do not vanish. These are as follows(see middle panel of Fig. 9):

• 1.

  The dark vertical stripe at g − r = 0.2 is a sign of too many blue thin-disc stars (A stars), which do not show up in Figs. 8 and 7. The constant SFR overestimates the fraction of young stars.

• 2.

  The white island at g − r = 0.35 and 15 < g < 17 mag points to a mismatch of the thick disc turn-off stars. There are too few nearby stars owing to the flat sech² profile and the turn-off is too red.

• 3.

  The faint red end with a colour–magnitude range of g − r > 1.0 and g > 18 mag shows too many M dwarfs. The reason may be that the g − r colour of the thin disc M dwarfs is slightly too blue or that the IMF allows too many stars at the low mass end.

• 4.

  The white islands at the bright red end are another hint that the red giants of the thick disc and halo are too blue (see also Figs. 3 and 11).

The magnitude distribution in the g and r bands (see Fig. 10) shows larger discrepancies than for the default Set 1 (top row of Fig. 4). The shoulder at 17   mag < g < 19   mag is still missing. This is again a hint that the Hess diagram analysis is a much more powerful constraint of the stellar population than luminosity function fitting.

To quantify the impact of giants on the Hess diagram, we select the giants using the surface gravity log g < 3.8. The absolute and relative contributions of giants simulated by parameter set 5 are shown in Fig. 11. The giants of  ~30% stars exist in the Hess diagram. The rigid line of giant counts is close to g − r = 0.5 ~ 0.7. The giants in the thin disc have a higher relative contribution, though the absolute counts are too low for the thin disc.

The parameters of the optimized model are listed in Table 3 (set 5). Compared to the default model (Set 1), the final Set 5 has a lower local surface density, a smaller maximum scale height h_z,max of the thin disc, and a shallower nearby halo profile of lower local density, because the old, extended thin disc and the halo do not need to account also for the thick disc stars.

Compared to the Just-Jahreiß model with  ⟨ χ² ⟩  = 4.31 (Just et al. 2011), the fit of Set 5 is much worse. This is also obvious from the median relative deviations Δ_med = 0.2576. The mismatch of the isochrones in the ugriz filters alone cannot account for these deviations. The simple structural forms of the thin and thick discs are the main reasons for the large discrepancies.

5. Simulations and results of the Besançon model

The Besançon model allows no variations in the choice of free parameters. We performed a multi-colour analysis of the predictions in different colours, as given in Table 5. We selected all stars in the NGP field from DR7 within the colour and apparent magnitude ranges listed in Table 5⁴ and with magnitude errors smaller than 0.2 mag in each filter. The selection of the field and objects is similar to that used in Just et al. (2011) and Sect. 4.

5.1. Isochrones

Since the thick disc and the halo are modelled in the Besançon model by simple stellar populations with a single age and a relatively narrow metallicity range (11 Gyr and  [Fe/H]  = −0.78 ± 0.30 for the thick disc, 14 Gyr and  [Fe/H]  = −1.78 ± 0.50 for the halo), we can compare the isochrones with the fiducial sequences of star clusters with corresponding properties (see Table 1). In An et al. (2009), a detailed analysis of isochrone fitting to fiducial sequences in the ugriz filters of globular and open clusters from An et al. (2008) is presented. The authors included a comparison with clusters observed with the CFHT-MegaCam by Clem et al. (2008). An et al. (2009) found that the clusters show a systematic offset relative to the Yale Rotating Evolutionary Code (YREC) and MARCS models owing to a zero-point problem. We use five old globular clusters to cover the metallicity range of the thick disc and halo.

We plot CMDs using the absolute magnitudes M_g, M_r, and M_i from −3 mag to 15 mag and compare the CMDs to the fiducial isochrones of the globular clusters presented in Table 1. There is a considerable offset of 0.3 mag in (u − g) between the fiducial isochrones and the simulated CMDs (top row of Fig. 12), which is a hint of a zero-point problem in the u band. This large offset explains why the Hess diagram of the u − g completely fail to the actual data (see below).

In Fig. 12, the left panels show the comparison for the thick disc, whereas the right-hand panels are for the halo. In the different CMDs, the fiducial isochrones of the globular clusters are overplotted: M 5 and M 71 are used to represent the thick disc; and M 5, NGC 4147 and M 92 are used for the halo (plotted in sequence from the turn-off point from the left to right in the right panels).

In the other CMDs of the thick disc, there is a smaller, but still significant, offset in the blue colour for the MS in the Besançon model compared to the fiducial isochrones. The metallicity of the cluster M 71 is closest to the metallicity of the thick disc adopted in Robin et al. (2003). The matches of M 5 with the lower metallicity are obviously better than for M 71, though there is still an offset (~0.2 mag in the u − g) between the fiducial isochrones of M 5 and the sequence of the thick disc of the Besançon model. The isochrones of the thick disc correspond to fiducial isochrones with a metallicity below −1.0.

For the halo, the offsets are similar but slightly larger. The isochrone with  [Fe/H]  = −2.28 is closest to the model but corresponds to the lower boundary of the metallicity range in the Besançon model.

5.2. Luminosity functions

In Fig. 13, the luminosity functions of the SDSS data (with error bars) and the simulated data are compared in all five filters. The full (black) step-function shows all stars of the Besançon model with contributions from the thin disc (dotted red lines), thick disc (dashed cyan lines), and halo (dot-dashed green lines). The model fits the data well over a wide magnitude range in all filters. A few systematic features can be seen:

• in the u, g, i, and z bands, the predicted number densities are slightly higher than in the actual data. At the brighter magnitudes (g < 15 mag), these differences are significant in all five filters. The thin disc dominates at the bright end where the differences between the predictions and the data show the largest deviations (0.2 dex to 0.5 dex). These stars correspond to F and early G stars (M_g = 4 to 5 mag) of the MS of the thin disc with distances above  ~500 pc (see Sect. 5.3). The reason for the discrepancy may be a higher mass density of the thin disc at large z or a higher fraction of early-type stars caused by the IMF or the SFH;

• only in the r band, the star counts are underestimated at r > 19th mag;

• the thin disc dominates only at the bright end, but contributes significantly up to the faint magnitude limit in the r, i, and z bands;

• the contribution of the thick disc dominates over wide magnitude ranges showing that the thick disc model parameters are very important in attempting to reproduce the luminosity functions;

• the halo contributes significantly only at the faint end of the u, g, and r band diagrams.

5.3. Hess diagrams

We calculate the Hess diagrams of the SDSS data and the Besançon model in the same way. The star counts are binned in each CMD in intervals of 0.05 mag in colour and 0.5 mag in magnitude and smoothed in steps of 0.01 mag. Numerical counts are given on a logarithmic scale covering a range of 1 to 100 per deg², 0.1 mag in colour, and 1 mag in luminosity. The relative differences between the model and the data (i.e., (data−model)/model) are smoothed only in luminosity and cover a linear range of −1.0 to  + 1.0. We investigate star counts in the four CMDs (u − g,g), (g − r,g), (r − i,r), and (i − z,i). To show and compare each stellar population in all the CMDs, the colour and magnitude ranges are different in the different CMDs (see Cols. 2 and 3 of Table 5).

In Figures 14–17, the results are shown for all CMDs. The upper panels show the Hess diagrams of the data, the model, and their relative differences (from left to right). In the lower panels, the contributions of the thin disc, the thick disc, and the halo (from left to right) as given by the Besançon model are shown with the same colour coding.

We first discuss the results in (g − r,g) for a direct comparison to the Just-Jahreiß (Just & Jahreiß 2010) and the TRILEGAL (Girardi et al. 2005) models. The top right panel in Fig. 15 shows that the relative differences for most regions are smaller than 50%. There are four regions in the CMD exceeding that level:

• the deep minimum at the lower right corner, where the star count predictions exceed the observations by more than a factor of three (compare the middle and left panels in the top row) is due to M dwarfs in the thin disc model (lower left panel). It seems that the (g − r) colour of the M dwarfs is too blue;

• the F turn-off regions of the thick disc (lower middle panel) and the halo (lower right panel) are also blue-shifted by  ~0.1 mag in (g − r), leading to the negative-positive features at g = 16 and g = 20 mag, respectively, in the difference plot;

• at the bright end (g = 14 to 15 mag), the star counts are determined by MS stars of the nearby thin disc, as well as bright F turn-off stars and giants in the thick disc. The thick-disc giants are responsible for the excess of stars on the blue side of the diagram. The maximum at the top right corner is a sign of missing thin-disc K dwarfs or red giants in the thick disc;

• the faint blue part at g > 19 mag and g − r = 0.2 to 0.4 mag is dominated by the halo. The density of the halo is strongly under-represented by the model. This cannot be seen in the luminosity function (see Fig. 13), since it is compensated for by the over-estimated number of thin disc M dwarfs at g − r < 1.2;

• the Besançon model confirms that the faint, very blue plume at g − r < 0.2 cannot be reproduced by stellar populations of the Milky Way and is most probably dominated by mis-identified extragalactic sources (Lupton et al. 2001) (see also Schultheis et al. 2006).

In the Just-Jahreiß model, the median relative deviation in the Hess diagram is Δ_med = 5.63%. This corresponds to a reduced χ² value, measured by the star counts in logarithmic scale, of 4.31. In the Besançon model, typical deviations are five times larger with a corresponding χ² = 408.15. The median relative deviations of the different colour diagrams are between 20% and 53%. In the CMD (Hess diagram), the star count distribution at a given colour is mainly determined by the geometrical structure of the components, whereas the distribution at a given magnitude is strongly influenced by the colours of the underlying stellar populations. Since all strong features in the relative differences between the data and model in Fig. 15 show a horizontal structure, the main reason for this disagreement with in the Besançon model are the mis-matching properties of the synthetic stellar populations (see Fig. 12). In the Just-Jahreiß model, the MS properties of the stellar components fit by construction in the solar neighbourhood.

In the TRILEGAL model (Sect. 4), the colour offsets in the synthetic stellar populations are much smaller and do not completely dominate the deviations of the data from the model. For features in the Hess diagram other than the ones just discussed, the absolute relative difference has a median value of 0.2576 corresponding to χ² = 104.19.

In case of the Besançon model, the overall differences between the data and model arising from deviations of the density profiles are large and cannot be quantified here.

In (u − g,g) (see Fig. 14) diagram the differences between the data and the Besançon model are much larger than in the (g − r,g) plot. The main reason is that the isochrones are even more blue-shifted by  ~0.3 mag in (u − g) (see Fig. 12). In addition, the features in the number density distributions of the Hess diagrams are much wider in colour in the observations than in the model.

In the (r − i,r) (see Fig. 16) diagram, the thin disc dominates for r − i > 0.8, while the thick disc dominates for 0.3 < r − i < 0.7. That the largest significant difference between the model and the observations is interestingly located at the interface of the populations of the thin and thick disc is caused by a gap at 0.75 < r − i < 0.85 and 20 < r < 20.5. Another feature is a pronounced gap at r − i ≈ 0.45 in the MS of the thin disc (lower left panel), which is much weaker in the other colours. These minima correspond to the relative number counts along the MS of the thin disc. At red colours (low-mass end), the star numbers are determined by the IMF, whereas in the blue (high-mass end) there is a strong impact of the SFH owing to the finite MS lifetime of the stars.

The Hess diagram of (i − z,i) (Fig. 17) represents the best fit of the Besançon model and SDSS data. In most regions of the CMD the median absolute difference is about 0.20.

Altogetherour multi-colour analysis confirms that the synthetic population properties first of all need to be improved in the Besançon model. In a second step, the IMF, SFH, and the spatial structure should be investigated with higher precision.

6. Summary

We have compared star-count predictions generated by the TRILEGAL and Besançon models with the observed sample in the NGP field (b > 80°) from the SDSS DR7 photometric stellar catalogue using de-reddened magnitudes adopting the extinction values of Schlegel et al. (1998). We measured the quality of the agreement between different models and the data by performing a χ² determination of the Hess diagrams, i.e., the star count distribution in the CMD. As a comparison, we referred to the Just-Jahreiß model (model A in Just & Jahreiß 2010; Just et al. 2011) with  ⟨ χ² ⟩  = 4.31 corresponding to a median relative deviation between the data and model of Δ_med = 5.63% in (g − r,g). Additionally, fiducial isochrones of five globular clusters, based on the SDSS observations of An et al. (2008) combined with the cluster parameters of Harris (1996), were used to test the synthetic stellar populations of the thick disc and halo modelled by the TRILEGAL and Besançon models.

6.1. The TRILEGAL model

We studied with TRILEGAL five different models simulating the NGP region and compared their output to observational photometric data of this region from the SDSS. The web interface of TRILEGAL allowed us a restricted choice of input parameters and functions for the Galactic components. We compared the default input set (lacking a thick disc component) with best-fit adjustments for three other Galaxy models from the literature, namely the Just-Jahreiß model (Just & Jahreiß 2010), the Besançon model (Robin et al. 2003), and the model of Jurić et al. (2008). We found that the default set reproduces best the luminosity functions in the g and r bands. The  ⟨ χ² ⟩  = 151.98 from the Hess diagram analysis is slightly worse than that based on the Jurić model with  ⟨ χ² ⟩  = 137.88. All these models show a typical discrepancy Δ_med of 26% to 46% in the Hess diagram, which is significantly larger than the disagreements from the Just-Jahreiß model.

A comparison of the colour–absolute magnitude diagrams of the thick disc and halo with the observed fiducial isochrones of globular clusters (An et al. 2008) shows that the MSs and the F-turnoff regions of the stellar evolution library are roughly consistent in the g and r filters. The small colour offset of g − r ~ 0.05 mag may be due to age differences among the components. However, the giant branches show large systematic disagreements. The red colour of the K giants observed in the globular clusters are not reproduced.

Since the parameter space was too large to find the absolute minimum  ⟨ χ² ⟩  with TRILEGAL, we determined the parameters for halo, thick disc, and thin disc sequentially by fitting smaller areas of the CMD. Finally we found an optimised set 5 with  ⟨ χ² ⟩  = 1.33 and a median deviation Δ_med of 26%, which is five times larger than in the Just-Jahreiß model. The parameters of this model are listed in the third-last column of Table 3. The reason for the larger deviations is mainly the restricted choice of thin and thick disc parameters (for a constant SFR and the adopted vertical density-profile shapes).

The TRILEGAL code is a powerful tool for testing Milky Way models against star count data. The investigation of star number densities in the CMD (Hess diagrams) are particularly helpful in distinguishing the impacts of the different stellar populations.

Some more progress may be possible by varying the thin disc parameters further, but the main disagreements are due to fundamental limitations of the models. In the TRILEGAL model, the range of input parameters and functions is quite restricted. More freedom in choosing the SFH, the vertical profiles of the thin and thick discs and for halo profile would be desirable when simulating high Galactic-latitude data. The comparison with fiducial star cluster isochrones has additionally shown that an improvement in the predicted colours of stars for late stages of stellar evolution in the ugriz filter system is desirable.

6.2. The Besançon model

We investigated in a multi-colour analysis the properties of the Besançon model by comparing with SDSS data of the NGP field. We generated mock stellar catalogues with the Besançon model. The output in u^∗, g′, r′, i′, and z′ photometries was converted into the SDSS u, g, r, i, and z system using the linear relations provided by Tucker et al. (2006) and Just & Jahreiß (2008).

The comparison of the colour–absolute magnitude CMDs of the model output of thick disc and halo with the fiducial globular clusters covering the respective metallicity ranges, uncovered significant differences. The MS and turnoff colours are blue-shifted by from 0.05 mag to 0.3 mag, increasing from the red to the blue colour bands. The large offset in u − g is mainly caused by a zero-point problem in the u′ filter of the CFHT-MegaCam reported by An et al. (2008). The giant branches are fitted well but the models do not extend far enough to the red end of the CMDs.

A comparison of the luminosity functions in the different filters shows that there is a good general agreement between the data and the model. Only for magnitudes  < 15 mag, the star counts are overestimated. In the r band, there is a deficit of stars at the faint end at r > 19 mag. The contributions of the different components (thin disc, thick disc, and halo) as a function of apparent magnitude were plotted to identify the dominating populations. The fraction of giants are generally of the order of a few percent, which may be neglected for the sake of our tests.

The Besançon model output and the differences from the data were plotted in Hess diagrams of four colours [(u − g,g), (g − r,g), (r − i,r), and (i − z,i)]. To analyse the impact of the different Galactic components, the individual Hess diagrams of the three components were also shown. The relative difference between the data and the model cover a range larger than −100% to 100% in each Hess diagram. Except for some well-defined features caused by the colour offset of the isochrones, the median relative differences in the Hess diagrams are between 0.2 and 0.5 (ignoring the strongly perturbed (u − g,g) diagram). The largest discrepancies in the Hess diagrams are the follows:

• the MS turn-off regions of thick disc and halo are underestimated in terms of the observed colours, meanwhile they are overestimated in terms of the modelled colours owing to the blue-shifted MS in the Besançon model;

• the halo is strongly underestimated, which is compensated for in the luminosity functions by M dwarfs of the thin disc, which fall in the model in the investigated colour range but do not do so in the data;

• the excess of observed stars in the transition between the thin disc and the thick disc compared to the model in terms of blue colours may also be due to the shifts in the synthetic colours produced by;

• the prominent gap in the MS of the thin disc at r ~ 0.45, which is not present in the SDSS data. This feature appears in all of the thin disc realizations of the Hess diagrams, but the gap in the r − i diagram is the most significant one.

6.3. General conclusions

We have investigated the ability of the TRILEGAL model to reproduce star counts in the SDSS filters (g − r,g) and the quality of the Besançon model in the full-filter system ugriz for the NGP field. We have analysed the number distributions of stars in the CMDs (via Hess diagrams) instead of only using luminosity functions. We have used the Just-Jahreiß model (Just & Jahreiß 2010; Just et al. 2011) to quantify the quality of the fits. In comparison to the Just-Jahreiß model with a median deviation between data and model of 5.6%, the TRILEGAL model with deviations of the order of 26% and the Besançon model with median deviations reaching from 20% to 53% show significant deficits. In the TRILEGAL model, the main reason for the larger deviations is the restricted choice of the density profiles for the thin disc, thick disc, and halo. In the Besançon model, the synthetic colours in the ugriz filters of the stellar populations need first of all to be improved, before the other properties of the Milky Way components can be investigated further.

In general, Milky Way models based on synthetic stellar populations such as the TRILEGAL and Besançon models are powerful tools for producing mock catalogues of stellar populations for comparison with data. Both models, TRILEGAL and Besançon, reproduce the observed SDSS luminosity functions very well. However, the luminosity functions are insensitive to the assumed structural properties of the Milky Way components. To improve the modelled Hess diagrams, the synthetic colours of the stellar populations in ugriz filters need to be upgraded in both the TRILEGAL and the Besançon models. We are looking

forward to applying the TRILEGAL, the Besançon, or the new Galaxia (Sharma et al. 2011) tool to the full SDSS data set.

---

Acknowledgments

We thank the anonymous referee for thoroughly reading the paper and providing comments that helped us to improve the manuscript. We are grateful to Dr. L. Girardi for important discussions about the TRILEGAL model, and to Dr. S. Vidrih for his contributions to the codes used for the comparisons presented in this paper. S.G. thanks Dr. A. Robin for her discussions on the local calibrations of the Besançon model. S.G. was supported by a scholarship of the China Scholarship Council (CSC). This research was also supported by the Collaborative Research Centre (Sonderforschungsbereich SFB) 881 “The Milky Way System”. Funding for the SDSS and SDSS-II has been provided by the Alfred P. Sloan Foundation, the Participating Institutions, the National Science Foundation, the US Department of Energy, the National Aeronautics and Space Administration, the Japanese Monbukagakusho, the Max Planck Society, and the Higher Education Funding Council for England. The SDSS web site is http://www.sdss.org/. The SDSS is managed by the Astrophysical Research Consortium for the Participating Institutions. The Participating Institutions are the American Museum of Natural History, Astrophysical Institute Potsdam, University of Basel, University of Cambridge, Case Western Reserve University, University of Chicago, Drexel University, Fermilab, the Institute for Advanced Study, the Japan Participation Group, Johns Hopkins University, the Joint Institute for Nuclear Astrophysics, the Kavli Institute for Particle Astrophysics and Cosmology, the Korean Scientist Group, the Chinese Academy of Sciences (LAMOST), Los Alamos National Laboratory, the Max-Planck-Institute for Astronomy (MPIA), the Max-Planck-Institute for Astrophysics (MPA), New Mexico State University, Ohio State University, University of Pittsburgh, University of Portsmouth, Princeton University, the United States Naval Observatory, and the University of Washington.

References

All Tables

All Figures