The two-dimensional displacement joint probability density P-Delta(X,Z) for water flowing through a bed of glass beads has been measured by means of pulsed field-gradient nuclear magnetic resonance. The simultaneous particle displacements X and Z perpendicular and parallel to the pressure gradient, respectively, at a given encoding time Delta, are obtained from an experiment employing orthogonal magnetic field gradients. The resulting probability density distribution is compared to numerical simulations of flow through an equivalent system of randomly deposited monodisperse spheres. The dependence of the centered second moments in X and Z on flow time is discussed for the experimental and simulated data. A crossover from a time scale dominated by Brownian motion toward a behavior determined by the convective flow and velocity fluctuations is observed. The mutual dependence between displacements perpendicular and parallel to the flow direction is revealed in the evolution of a correlation coefficient rho(X2,Z) This coefficient is found to increase for short times and to decrease for larger displacements, with a maximum at an average displacement corresponding to the bead radius. As a means of displaying the cause of these correlations, a correlation probability density C-Delta(X,Z)= P-Delta(X,Z) - P-Delta(X)P-Delta(Z) is suggested, where P-Delta(X) and P-Delta(Z) are the marginals of P-Delta(X,Z). A plot of this matrix renders zero in the absence of correlations, but produces a characteristic pattern of positive and negative regions when displacements in X are correlated with those in Z. The time evolution of this pattern is discussed and compared to the shape of model propagators obtained from an analytical function and a numerical simulation for a simplified capillary array, respectively. [S1063-651X(98)05611-6].