Natural neighbor interpolation is a weighted average interpolation method that is based on Voronoi tessellation. In this paper, we present and implement an algorithm for performing natural neighbor interpolation using graphics hardware. Unlike traditional software-based approaches that process one query at a time, we develop a scheme that computes the entire scalar field induced by natural neighbor interpolation, at which point a query is a trivial array lookup, arid range queries over the field are easy to perform. Our approach is faster than the best known software implementations and makes use of general purpose stream programming capabilities of current graphics cards. We also present a simple scheme that requires no advanced graphics capabilities and can process natural neighbor queries faster than existing software-based approaches. Finally, recognizing the limitation incurred by the bounded size of graphics frame buffers, we propose a sub-division approach that allows performing queries locally in a subdivision of the input domain. This approach can reduce to a negligibly small degree (< 1%) the loss of precision caused by the naive scaling method while still processing queries faster than the software-based approaches when the number of sites is large.