Terrestrial carbon processes are affected by soil moisture, soil temperature, nitrogen availability and solar radiation, among other factors. Most of the current ecosystem biogeochemistry models represent one point in space, and have limited characterization of hydrologic processes. Therefore, these models can neither resolve the topographically driven spatial variability of water, energy, and nutrient, nor their effects on carbon processes. A spatially-distributed land surface hydrologic biogeochemistry model, Flux-PIHM-BGC, is developed by coupling the Biome-BGC model with a physically-based land surface hydrologic model, Flux-PIHM. In the coupled system, each Flux-PIHM model grid couples a 1-D Biome-BGC model. In addition, a topographic solar radiation module and an advection-driven nitrogen transport module are added to represent the impact of topography on nutrient transport and solar energy distribution. Because Flux-PIHM is able to simulate lateral groundwater flow and represent the land surface heterogeneities caused by topography, Flux-PIHM-BGC is capable of simulating the complex interaction among water, energy, nutrient, and carbon in time and space. The Flux-PIHM-BGC model is tested at the Susquehanna/Shale Hills Critical Zone Observatory. Model results show that distributions of carbon and nitrogen stocks and fluxes are strongly affected by topography and landscape position, and tree growth is nitrogen limited. The predicted aboveground and soil carbon distributions generally agree with the macro patterns observed. Although the model underestimates the spatial variation, the predicted watershed average values are close to the observations. The coupled Flux-PIHM-BGC model provides an important tool to study spatial variations in terrestrial carbon and nitrogen processes and their interactions with environmental factors, and to predict the spatial structure of the responses of ecosystems to climate change.