Forming a symmetric, unreduced, tridiagonal matrix with a given spectrum

Given a set of $n$ distinct real numbers, our goal is to form a symmetric, unreduced, tridiagonal, matrix with those numbers as eigenvalues. We give an algorithm which is a stable implementation of a naive algorithm forming the characteristic polynomial and then using a technique of Schmeisser.