In this paper, we present a detailed next-to-leading-order (NLO) study of $J/psi$ angular distributions in $e^{+}e^{-} to J/psi+eta_c,chi_{cJ}$ ($J=0,1,2$) within the nonrelativistic QCD factorization (NRQCD). The numerical NLO expressions for total and differential cross sections, i.e., $frac{dsigma}{dcostheta}=A+Bcos^2theta$, are both derived. With the inclusion of the newly-calculated QCD corrections to $A$ and $B$, the $alpha_{theta}(= B/A)$ parameters in $J/psi+chi_{c0}$ and $J/psi+chi_{c1}$ are moderately enhanced, while the magnitude of ${alpha_theta}_{J/psi+chi_{c2}}$ is significantly reduced; regarding the production of $J/psi+eta_c$, the $alpha_theta$ value remains unchanged. By comparing with experiment, we find the predicted ${alpha_theta}_{J/psi+eta_c}$ is in good agreement with the $textrm{B}scriptsize{textrm{ELLE}}$ measurement; however, ${alpha_theta}_{J/psi+chi_{c0}}$ is still totally incompatible with the experimental result, and this discrepancy seems to hardly be cured by proper choices of the charm-quark mass, the renormalization scale, and the NRQCD matrix elements.