A digital surface model (DSM) is an important geospatial infrastructure used in various fields. In this paper, we deal with how to improve the quality of DSMs generated from stereo image matching. During stereo image matching, there are outliers due to mismatches, and non-matching regions due to match failure. Such outliers and non-matching regions have to be corrected accurately and efficiently for high-quality DSM generation. This process has been performed by applying a local distribution model, such as inverse distance weight (IDW), or by forming a triangulated irregular network (TIN). However, if the area of non-matching regions is large, it is not trivial to interpolate elevation values using neighboring cells. In this study, we proposed a new DSM interpolation method using a 3D mesh model, which is more robust to outliers and large holes. We compared mesh-based DSM with IDW-based DSM and analyzed the characteristics of each. The accuracy of the mesh-based DSM was a 2.80 m root mean square error (RMSE), while that for the IDW-based DSM was 3.22 m. While the mesh-based DSM successfully removed empty grid cells and outliers, the IDW-based DSM had sharper object boundaries. Because of the nature of surface reconstruction, object boundaries appeared smoother on the mesh-based DSM. We further propose a method of integrating the two DSMs. The integrated DSM maintains the sharpness of object boundaries without significant accuracy degradation. The contribution of this paper is the use of 3D mesh models (which have mainly been used for 3D visualization) for efficient removal of outliers and non-matching regions without a priori knowledge of surface types.