Given scattered data in IR s , interpolation from a dilated box spline space S M (2 k ·)is always possible for a fine enough scaling. For example, for the Lagrange function of a point θ one could take any shifted dilate M (2 k · −j) which is nonzero at θ and zero at the other interpolation points. However, the resulting interpolant, though smooth (and local), will consist of a set of "bumps", and so by any reasonable measure provides a poor representation of the shape of the underlying function. On the other hand, it is possible to choose a space of interpolants which contains some M (2 k · −j) of arbitrarily large support. But the resulting methods are increasingly less local, and in general still require some splines with a much higher level of dilation. Here we provide a multilevel method which constructs a space of interpolants by taking as many splines as possible from a given dilation level, then as many from the next (higher) dilation level, and so forth. The choice at each level is made using the suggestion of [W99], which is based on the Riesz representation theorem. This requires an inner product on the ground space S M , and the higher levels. .. The inner products used here involve the box spline coefficients, and prewavelet coefficients of [RS92], respectively, and are norm equivalent to · L 2 (IR s ) . These lead to a scheme which is easily implemented, and numerically stable. Previously, box spline interpolants have been considered only for data on a regular grid.