Bulk-density ( ρ) of soil is an important indicator of soil compaction and type. A knowledge of the spatial variability of in situ soil density is important in geotechnical engineering, hydrology and agriculture. Surface geophysical methods have so far shown limited success in providing an accurate and high-resolution image of 3-D soil-density distribution. In this pursuit, 3-D seismic full-waveform inversion (FWI) is promising, provided the robustness and accuracy of density inversion via this approach can be established in the near-surface scale. Ho wever , simultaneous reconstruction of ρand seismic wave velocities through multiparameter FWI remains a challenging task. Near-surface seismic data are commonly dominated by dispersive surface waves whose velocities are controlled by the value and distribution of shear-wave velocity ( V S ). One major difficulty in estimating reliab ly ρfrom near -surface seismic data is due to the relati vel y low sensiti vity of the seismic w av efield to ρcompared to seismic v elocities. Additionally, the accuracy of the estimated ρdecreases due to error in V S -an issue known as parameter coupling. Parameter coupling makes it difficult to estimate accurately ρwithin the framework of conventional gradient-based FWI. More sophisticated optimization approaches (e.g. truncated Newton) can reduce the effect of parameter coupling, but these approaches are commonl y not af fordab le in near -surface applications due to heavy computational burden. In this research, w e ha v e inv estigated how choosing correctly the force direction of the seismic source can contribute to a higher accuracy of ρestimates through 3-D FWI. Using scattered wavefields, the Hessian, and inversion tests, an in-depth and systematic investigation of data sets corresponding to different force directions has been carried out. A comparison of the scattered wavefields due to a point-localized ρperturbation for different force directions shows the robustness of the horizontal-force data set to noise compared to the vertical-force data set. Fur ther more, for a point-scatterer model, an analysis of the gradients of the misfit function using the Hessian shows that utilizing a horizontal-force source enables one to reconstruct the high-resolution gradient with relati vel y small parameter coupling. Finally, inversion tests for two different subsoil models demonstrate that 3-D FWI on a horizontal-force-source seismic data set is capable of providing a more accurate 3-D ρdistribution in soil compared to a vertical-force-source data set. Our results show that the use of a horizontal-force source might allow avoiding computationally demanding, costly optimization approaches in 3-D FWI.