We implement the universal wave-function overlap (UWFO) method to extract modular
S and
T matrices for topological orders in Gutzwiller-projected parton wave functions (GPWFs). The modular
S and
T matrices generate a projective representation of
SL(2,Z) on the degenerate-ground-state Hilbert space on a torus and may fully characterize the 2+1D topological orders, i.e., the quasiparticle statistics and chiral central charge (up to
E8 bosonic quantum Hall states). We use the variational Monte Carlo method to computed the
S and
T matrices of the chiral spin liquid (CSL) constructed by the GPWF on the square lattice, and we confirm that the CSL carries the same topological order as the $\ensuremath{\nu}=\frac{1}{2}$ bosonic Laughlin state. We find that the nonuniversal exponents in the UWFO can be small, and direct numerical computation can be applied on relatively large systems. The UWFO may be a powerful method to calculate the topological order in GPWFs.