We present a generalized and statistically coherent scheme of computing hardness ratios and associated error bars. In this scheme, we model the observed counts as a non-homogeneous Poisson process and exploit sophisticated Bayesian approaches (e.g., Gibbs sampling) to calculate hardness ratios, accounting for local background contamination and effective area variations. We apply this scheme to the simple counts ratio [S/H] as well as its variants, colors [log(S/H)] and fractional difference hardness ratios [(H-S)/(H+S)]. We also perform simulations to compare the new Bayesian methods with the classical method, thereby illustrating that (a) the former provides more accurate estimates of the uncertainties, (b) the mode of the posterior probability distribution function (pdf) is a robust estimator of the hardness ratio, and (c) the pdfs of the colors are the best behaved in the low counts limit.

We apply this method to identify candidate qLMXB's in the globular cluster Terzan 5.