Abstract The motion of baryonic components of the Milky Way is governed by both luminous and dark matter content of the Galaxy. Thus, the dynamics of Milky Way globular clusters (GCs) can be used as tracers to infer the mass model of the Galaxy up to a large radius. In this work, we apply the directly observable line-of-sight velocities to test if the dynamics of the GC population are consistent with an assumed axisymmetric gravitational potential of the Milky Way. For this, we numerically compute the phase space distribution of the GC population where the orbits are either oriented randomly or co-/counter- rotating with respect to the stellar disk. Then we compare the observed position and line-of-sight velocity distribution of ∼150 GCs with those of the models. We found that, for the adopted mass model, the co-rotating scenario is the favored model based on various statistical tests. We do the analysis with and without the GCs associated with the progenitors of early merger events. This analysis can be extended in the near future to include precise and copious data to better constrain the Galactic potential up to a large radius.