Grain size distribution (GSD) is crucial for understanding soil properties and surface processes. We find that both terrestrial soils and lunar soils are subjected to a unified GSD function, P(D)= g(mu )D-mu exp(-D/Dc), reducing the textural fractions and grade modes to a parameter pair (mu, Dc), which unifies terrestrial and extraterrestrial soils in granular configuration, beyond the environments and mechanisms of soil genesis. To construct a framework of the soil formation, we generalize the textural composition to a grade space representing the granular configuration, and conceptualize soil genesis as the random aggregation of the fractal fragmentation of parent lithospheric material and fragments from other sources (e.g., meteorites impacts or surface transport processes). Random simulation reproduces the multiple grade modes observed in soils, and spontaneously derives the unified GSD function. Then we numerically generate the (mu, Dc)-fields for soils on earth and moon, which refine the digital data mapping based on site measurements and depict the local fluctuation of soil parameters. The GSD unity also provides a tool of generating "numerical simulants" of lunar soils to fill the gap in material simulants. The study leads to a GSD-paradigm (in contrast to the conventional landscape-paradigm) in soil study, which is expected to facilitate the data harmonization on earth and promote the generation of lunar regolith data in favor of the in-situ resource utilization and base construction on moon.