This research develops a novel framework for interpolation and approximation of multivariate data. The research introduces a method for construction of non-separable multivariate splines that are geometrically tailored for general sampling lattices with emphasis on sphere packing and covering lattices in multidimensions. These multivariate splines lead to the design of multidimensional signal processing tools, that have proven degree of continuity and approximation order for sampling lattices. The splines under investigation, called Voronoi splines, are B-spline-like elements that inherit the geometry of any sampling lattice from its Voronoi cell and form basis for reconstruction. Furthermore, the relationship of the proposed splines with the well-known multivariate box splines are investigated and used for evaluation purposes. Compared to other polyhedral splines, Voronoi splines form shiftinvariant spaces for approximation. This property makes them uniquely suitable for multidimensional signal processing.