We describe a new algorithm for scattered data interpolation. It is based on a modified Shepard method similar to that of Algorithm 660 but uses 10-parameter cosine series nodal functions in place of quadratic polynomials. Also, the interpolant has continuous second partial derivatives. An accompanying survey article presents test results that show the method to be more accurate than polynomial-based methods in terms of reproducing test functions with large variations and steep gradients.