Context. Regular follow-up of imaged companions to main-sequence stars often allows a projected orbital motion to be detected. Markov chain Monte Carlo (MCMC) has become very popular recent years for fitting and constraining their orbits. Some of these imaged companions appear to move on very eccentric, possibly unbound orbits. This is, in particular, the case for the exoplanet Fomalhaut b and the brown dwarf companion PZ Tel B on which we focus here. Aims. For these orbits, standard MCMC codes that assume only bound orbits may be inappropriate. Our goal is to develop a new MCMC implementation that is able to handle both bound and unbound orbits in a continuous manner, and to apply this to the cases of Fomalhaut b and PZ Tel B. Methods. We present here this code, based on the use of universal Keplerian variables and Stumpff functions. We present two versions of this code, the second one using a different set of angular variables that were designed to avoid degeneracies arising when the projected orbital motion is quasi-radial, as is the case for PZ Tel B. We also present additional observations of PZ Tel B. Results. The code is applied to Fomalhaut b and PZ Tel B. We confirm previous results in relation to, but we show that on the sole basis of the astrometric data, open orbital solutions are also possible. The eccentricity distribution nevertheless still peaks around ∼0.9 in the bound regime. We present a first successful orbital fit of PZ Tel B, which shows in particular that, while both bound and unbound orbital solutions are equally possible, the eccentricity distribution presents a sharp peak very close to e = 1, meaning a quasi-parabolic orbit. Conclusions. It has recently been suggested that the presence of unseen inner companions to imaged ones may lead orbital fitting algorithms to artificially give very high eccentricities. We show that this caveat is unlikely to apply to Fomalhaut b. Concerning PZ Tel B, we derive a possible solution, which involves an inner ∼12 M Jup companion, that would mimic a e = 1 orbit, despite a real eccentricity of around 0.7, but a dynamical analysis reveals that this type of system would not be stable. We thus conclude that our orbital fit is robust.