[Aims], Forest vegetation variability may be explained by the complex interplay among several spatial structuring factors, including climate and topography. We modelled the spatial variability of forest vegetation assemblages and significant environmental variables along a complex environmental gradient or coenocline to produce a detailed cartographic database portraying the distribution of forests along it. [Methods], We combined an analysis of ordination coenoclines with kriging over 772 field data plots from the third Spanish National Forest Inventory in an Atlantic–Mediterranean transitional area (northern Spain). [Important Findings], The best fitted empirical semivariogram revealed a strong spatial structure of forest species composition along the complex environmental gradient considered (the climatic–topographic gradient from north to south). The steady and gradual increase of semivariance with a marked lag distance indicates a gradual turnover of forest assemblages according to the climatic–topographic variations (regional or local). Two changes in the slope of the semivariogram suggest the existence of two different scales of spatial variation. The interpolation map by Kriging of forest vegetation assemblages along the main coenocline shows a clear spatial distribution pattern of trees and shrubs in accordance with the spatial variation of significant environmental variables. We concluded that the multivariate geostatistical approach is a suitable technique for spatial analysis of forest systems employing data from national forest inventories based on a regular network of field plots. The development of an assortment of maps describing changes in vegetation assemblages and variation in environmental variables is expected to be a suitable tool for an integrated forest management and planning.