A procedure is presented for the possible systematic development of exchange-correlation functionals using ab initio electron densities and accurate total energies. For a training set of first row open-and closed-shell systems, densities are computed and are used to determine asymptotically vanishing exchange-correlation potentials. The new functional is then written as an expansion in products of the density and its gradient, and optimum expansion parameters are determined through a least squares fit involving both these potentials and accurate exchange-correlation energies. Unlike conventional functionals, the potential of the fitted functional approaches a non-zero value asymptotically, and this is achieved by introducing a self-consistently computed system-dependent shift into the fitting procedure. This shift represents the influence of the integer derivative discontinuity in the exact energy. The method has been used to determine a 21 term spin-polarized exchange-correlation functional using Brueckner Doubles or MP2 densities of 20 small systems. For those with open-shells the computed shifts are close to the hardness of the system, while for closed-shells they are considerably smaller than the hardness. These observations are consistent with theoretical requirements. A comparison of the new potential with conventional potentials highlights important differences in the inter-shell and asymptotic regions, while the values of the shifts and highest occupied self-consistent eigenvalues suggest improved asymptotic densities. The mean absolute errors in self-consistent total energies and optimized bond-lengths of systems in the training set are 0.003E(h) and 0.01 Angstrom, respectively. Comparable values are obtained for 12 first-row closed-shell systems outside the training set. Compared to conventional functionals, the new functional predicts a significantly improved classical barrier height for the hydrogen abstraction reaction H+H-2-->H-2+H.