We investigate the low transverse momentum heavy quarkonium (J/ψ, Υ) productions in the forward rapidity region of pp and pA collisions at the LHC as an important probe to the transverse momentum tomography of the gluons in hadrons in the Color Glass Condensate (CGC) framework. By implementing the Sudakov resummation consistently in the CGC formalism, we achieve an excellent agreement between the improved CGC calculations and the LHC data. We show that both the small-x and the Sudakov effects are essential for a complete description of heavy quarkonium productions in the low transverse momentum region. This provides a solid foundation to study the small-x gluon saturation in a big nucleus with the future pA programs at the LHC.PACS numbers: 24.85.+p, 14.40.Pq, 12.38.Bx Introduction-Due to the rich and complex dynamics of strong interaction and its close connection to various interesting phenomena, heavy quarkonium production in high energy collisions has been one of the most intriguing subjects in the study of strong interaction physics. In particular, as the objective of this paper, it has been shown that the low-p ⊥ spectrum of forward rapidity heavy quarkoniums produced in high energy pp (pA) collisions carries important information of the gluon saturation at small-x. The concept of gluon saturation [1,2] comes from the anticipation of gluon recombination when its density becomes extremely high in high energy collisions. The associated dynamics is one of the challenging physics topics in current and future nuclear science program around the world [3].There have been intensive theoretical studies on heavy quarkonium productions in the small-x saturation framework, i.e., the Color-Glass-Condensate (CGC) formalism. They have shown that the typical transverse momentum of produced quarkonium is around the saturation momentum Q s (x g ) [4][5][6][7][8][9][10][11][12] . As a function of the longitudinal momentum fraction of target hadron x g , the saturation momentum increases when x g decreases, and it characterizes the typical transverse momentum that small-x gluons carry in a hadron. It can also be derived from the boundary which separates the dilute regime from the dense saturated regime where gluons evolve non-linearly. The non-linear evolution equation at small-x (or high energy) is known as the Balitsky-Kovchegov (BK) equation [13,14], which has been implemented numerically in CGC formalism.The physical picture of quarkonium productions in the CGC framework is as follows. Either before or after interacting with the dense gluons in the target hadron, a large x gluon from projectile proton fluctuates into a pair of heavy quarks, which eventually forms a bound state of quarkonium with transverse momentum p ⊥ . In the low p ⊥ region, the main contribution of the transverse momentum broadening is due to gluon saturation effects. At forward rapidity y, the kinematics indicates that x g from the target side is much smaller than the one in the projectile. This is the ideal region for the study of small-x physics. ...