This paper is concerned with the problem of computing the image of a set by a polynomial function. Such image computations constitute a crucial component in typical tools for set-based analysis of hybrid systems and embedded software with polynomial dynamics, which found applications in various engineering domains. One typical example is the computation of all states reachable from a given set in one step by a continuous dynamics described by a differential or difference equation. We propose a new algorithm for over-approximating such images based on the Bernstein representation of polynomial functions. The images are stored using template polyhedra. Using a prototype implementation, the performance of the algorithm was demonstrated on two practical systems as well as a number of randomly generated examples.