Calculate the state vector of planets and moons

What is the efficient way to calculate the state vector of planets and moons?

The coordinate system for planets would be with center in Solar System Baricenter, for Moon- with center in Earth and in Solar System Baricenter.

The use of JPL kernels (de430.bsp) and python jplephem package worked for me.
from jplephem.spk import SPK