The Levi-Civita wavefront coupling problem consists of, given two prescribed wavefronts, obtaining a refractive or reflective surface coupling them. We propose a robust numerical method to solve Levi-Civita’s problem, whose rationale is to consider that Levi-Civita’s solutions are level surfaces of the point characteristic function established between points of incoming and outgoing wavefronts. The method obtains both surface data points and their normals, enabling a more robust surface reconstruction. We carry out a detailed error analysis of our method by means of comparing the surface data estimation with nominal surfaces obtained in reference tests offering analytical solutions to Levi–Civita’s problem. The method offers, in computer simulations, highly accurate results with moderate computational cost.