The paper describes a fast algorithm for scattered data interpolation and approximation. Multilevel B-splines are introduced to compute a C/sup 2/ continuous surface through a set of irregularly spaced points. The algorithm makes use of a coarse to fine hierarchy of control lattices to generate a sequence of bicubic B-spline functions whose sum approaches the desired interpolation function. Large performance gains are realized by using B-spline refinement to reduce the sum of these functions into one equivalent B-spline function. Experimental results demonstrate that high fidelity reconstruction is possible from a selected set of sparse and irregular samples.