We present a new shape reconstruction method from photometric stereo. The method gives more accurate shapes than the existing method that uses three images and an algorithm based on the Frankot-Chellappa's to convert the surface normal to height. The new method, which is applicable to the case of two or more images available, solves the combination of two or more nonlinear equations representing the consistency between the reflectance maps and the observed images in an iterative way for the shape. And this is repeated for each of the four approximations of the surface normal and the results are averaged, to get a better shape estimate. Computer experiments are made to show the usefulness of the method and to present some of the characteristics.