Motivated by a recent halftoning method which is based on electrostatic principles, we analyze a halftoning framework where one minimizes a functional consisting of the difference of two convex functions (DC). One of them describes attracting forces caused by the image gray values, the other one enforces repulsion between points. In one dimension, the minimizers of our functional can be computed analytically and have the following desired properties: the points are pairwise distinct, lie within the image frame and can be placed at grid points. In the two-dimensional setting, we prove some useful properties of our functional like its coercivity and propose to compute a minimizer by a forwardbackward splitting algorithm. We suggest to compute the special sums occurring in each iteration step by a fast summation technique based on the fast Fourier transform at non-equispaced knots which requires only O(m log m) arithmetic operations for m points. Finally, we present numerical results showing the excellent performance of our dithering method.