Based on the nonrelativistic QCD factorization scheme, we present the first complete next-to-leadingorder study on the yield and polarization of Çð1S; 2S; 3SÞ hadroproduction. By using the color-octet longdistance matrix elements obtained from fits of the experimental measurements on Ç yield and polarization at the Tevatron and LHC, our results can explain the measurements on the yield very well, and for the polarizations of Çð1S; 2S; 3SÞ, they are in (good, good, bad) agreement with recent CMS measurement, but still have some distance from the CDF measurement.The well-known J=c polarization puzzle became obvious when the CDF measurement at the Tevatron [1] was found completely different from the leading-order (LO) theoretical prediction in the framework of nonrelativistic QCD (NRQCD) [2], which was proposed as a factorization approach on heavy quarkonium decay and production [3]. Even with the progress in the next-toleading-order (NLO) QCD calculation, theoretical studies [4][5][6] on J=c polarization at NLO could not clearly clarify the situation. Early measurements [7,8] on the polarization of Ç at the Tevatron are also in conflict with the corresponding LO NRQCD prediction [9]. Recently, a very important and interesting measurement on the polarization of Çð1S; 2S; 3SÞ at the LHC was reported by the CMS Collaboration [10], which employed improved consideration in the measurement [11]. Since the bottom quark is almost three times as heavy as the charm quark, the NRQCD s and velocity expansions are of better convergence and the theoretical predictions at QCD NLO are more reliable for Ç than that for J=c . Therefore, it is very important to extend the theoretical predication on Ç at QCD NLO to solve or clarify the long-standing polarization puzzle when there are already measurements at the LHC.In the last six years, there has been some very important progress in the NLO QCD correction calculation. The NLO corrections to color-singlet (CS) J=c hadroproduction have been investigated in Refs. [12,13], where its transverse momentum (p t ) distribution is found to be enhanced by 2-3 orders of magnitude at the high p t region and its polarization changes from transverse to longitudinal at NLO [13]. The results are reproduced at p t LO in a new factorization scheme for large p t quarkonium production [14]. The NLO corrections to J=c production via S-wave color-octet (CO) states ( 1 S ½8 0 , 3 S ½8 1 ) are studied in Ref.[15] and the corrections to the p t distributions of both the J=c yield and polarization are small. In Ref.[16], NLO corrections for cJ hadroproduction are studied. The complete NLO calculation for prompt J=c hadroproduction (with 3 P ½8 J included) was given by two groups [17,18], and their predictions for p t distributions agree with the experimental measurements at the Tevatron and LHC. The calculation for polarization of direct J=c hadroproduction at NLO QCD was presented by two groups [4,5]. The complete NLO calculation of the polarization for prompt J=c hadroproduction was completed by our group [6...