An extremely fast method is presented to calculate the local current–density vector in a flat conductor from the z-component of the magnetic field measured above its surface, e.g. by means of magneto-optical indicators, Hall-probe arrays or scanning SQUIDs. The method may be used for samples of arbitrary thickness provided that the current vector has only x- and y-components. The method combines the conjugate gradient (CG) method and fast Fourier transform to invert the relevant Toeplitz matrix equation. For a current map of n×n pixels, the number of operations needed is of order n2.8 only, compared to n4.5 or higher for earlier methods. The increase in speed for 512×512 pixels is found to be a factor 135 with respect to the fastest existing CG method.