  We present a new method based on a Bayesian hierarchical model to extract
constraints on cosmological parameters from SNIa data obtained with the SALT-II
lightcurve fitter. We demonstrate with simulated data sets that our method
delivers tighter statistical constraints on the cosmological parameters over
90% of the time, that it reduces statistical bias typically by a factor ~ 2-3
and that it has better coverage properties than the usual chi-squared approach.
As a further benefit, a full posterior probability distribution for the
dispersion of the intrinsic magnitude of SNe is obtained. We apply this method
to recent SNIa data, and by combining them with CMB and BAO data we obtain
Omega_m=0.28 +/- 0.02, Omega_Lambda=0.73 +/- 0.01 (assuming w=-1) and
Omega_m=0.28 +/- 0.01, w=-0.90 +/- 0.05 (assuming flatness; statistical
uncertainties only). We constrain the intrinsic dispersion of the B-band
magnitude of the SNIa population, obtaining sigma_mu^int = 0.13 +/- 0.01 [mag].
Applications to systematic uncertainties will be discussed in a forthcoming
paper.
