This paper deals with the structural analysis of rock fracture roughness, and accordingly, a method is developed for estimating/predicting the post-shearing 3D geometry of the fracture surface. For this purpose, surfaces of three natural rock fractures were digitized and studied before and after the direct shear test. The variogram analysis of the surfaces indicated a strong non-linear trend in the topography data. Hence, the spatial variability of the rock fracture surfaces was decomposed to: one deterministic component, characterized by a high-order polynomial representing the large-scale undulations, and one stochastic component, described by the variogram of residuals representing the small-scale roughness. Using an image-processing technique, a total of 343 damage zones with different sizes, shapes, initial roughness characteristics, local stress fields, and/or asperity strength values were spatially located and clustered. In order to characterize the overall spatial structure of the degraded zones, the concept of the ‘pseudo-zonal variogram’ was introduced. The results showed that the spatial continuity at the damage zones increases due to the asperity degradation. The increase in the variogram range is anisotropic and tends to be higher along the shearing. Consequently, the direction of maximum continuity rotates towards the shear direction. After modeling the evolution of the spatial structure with shearing and detecting boundaries of the degraded areas, a methodology was presented to provide a regression-kriging estimate of the morphology of sheared surfaces. The proposed method can be considered as a cost-free and reasonably accurate alternative to expensive techniques of scanning the rock fracture surface after the shear test.