We develop a mechanism to build the light-front wavefunctions (LFWFs) of meson bound states on a smallsized basis function representation. Unlike in a standard Hamiltonian formalism, the Hamiltonian in this method is implicit, and the information of the system is carried directly by the functional form and adjustable parameters of the LFWFs. In this work, we model the LFWFs for four charmonium states, η c , J/ψ, ψ , and ψ(3770) as superpositions of orthonormal basis functions. We choose the basis functions as eigenfunctions of an effective Hamiltonian, which has a longitudinal confining potential in addition to the transverse confining potential from light-front holographic QCD. We determine the basis function parameters and superposition coefficients by employing both guidance from the nonrelativistic description of the meson states and the experimental measurements of the meson decay widths. With the obtained wavefunctions, we study the features of those meson states, including charge radii and parton distribution functions. We use the J/ψ LFWF to calculate the meson production in diffractive deep inelastic scattering and ultra-peripheral heavy-ion collisions, and the η c LFWF to calculate its diphoton transition form factor. Both results show good agreement with experiments. The obtained LFWFs have simple-functional forms and can be readily used to predict additional experimental observables.