This article presents a method to produce accurate soil maps according to the numerical correlation of experimental data obtained from field description and soil laboratory analyses. At first, the diagnostic horizons and characteristics of all the 56 studied profiles were identified based on the standard definitions of the US soil classification system. In the next step, six thematic maps were prepared by interpolating the thickness of diagnostic horizons and/or the depth of the upper/lower boundaries of the soil layers with diagnostic characteristics. Afterward, the thematic maps were intersected to obtain the final soil map in this approach. The produced map with the present approach was compared with the conventional geopedologic soil map. The results showed that the geopedological approach could delineate 12 soil map units as consociation and 8 units as association while using the diagnostic layers as thematic maps were able to delineate all of the 20 soil map units as consociation. This approach introduced how the soil properties could be applied as a covariate in the SCORPAN model. Using interpolated soil, diagnostic layers could increase the precision and quality of soil maps.