  Following the operator definition of the fragmentation function developed by
Collins and Soper, we compute the gluon-to-$h_c$ fragmentation function at the
lowest order in the velocity expansion in NRQCD factorization approach.
Utilizing some modern technique developed in the area of multi-loop
calculation, we are able to analytically deduce the infrared-finite
color-singlet short-distance coefficient associated with the fragmentation
function. The fragmentation probability for gluon into $h_c$ is estimated to be
order $10^{-6}$.
