The lowest-lying glueball masses are computed in SU($N$) gauge theory on a spacetime lattice for constant value of the lattice spacing $a$ and for $N$ ranging from 3 to 8. The lattice spacing is fixed using the deconfinement temperature at temporal extension of the lattice $N_T = 6$. The calculation is conducted employing in each channel a variational ansatz performed on a large basis of operators that includes also torelon and (for the lightest states) scattering trial functions. This basis is constructed using an automatic algorithm that allows us to build operators of any size and shape in any irreducible representation of the cubic group. A good signal is extracted for the ground state and the first excitation in several symmetry channels. It is shown that all the observed states are well described by their large $N$ values, with modest $\cal O(1/N^2)$ corrections. In addition spurious states are identified that couple to torelon and scattering operators.