This paper presents a new method for the reconstruction of a surface from its ? and ? gradient field, measured, for example, via Photometric Stereo. The new algorithm produces the unique discrete surface whose gradients are equal to the measured gradients in the global vertical-distance least-squares sense. We show that it has been erroneously believed that this problem has been solved before via the solution of a Poisson equation. The numerical behaviour of the algorithm allows for reliable surface reconstruction on exceedingly large scales, e.g., full digital images; moreover, the algorithm is direct, i.e., non-iterative. We demonstrate the algorithm with synthetic data as well as real data obtained via photometric stereo. The algorithm does not exhibit a low-frequency bias and is not unrealistically constrained to arbitrary boundary conditions as in previous solutions. In fact, it is the first algorithm which can reconstruct a surface of polynomial degree two or higher exactly. It ...