This article presents a new approach to analytically compute high-frequency (HF) electromagnetic (EM) scattering from an arbitrary undulating rough surface according to the statistical characteristics of its geometric fluctuation. The integral expression of the scattering coefficient is derived from Kirchhoff approximation (KA) and extended to arbitrary spectra. Two methods are presented to calculate the integral: 1) the discrete Fourier transforming method based on the autocorrelation function (ACF) (KAC) and 2) the convolution method based on the power spectrum density function (KAW). For some analytically expressed spectra, the scattering coefficients can be calculated directly by the KAW and KAC methods according to its power spectrum density (PSD). The analytical expressions of classical KA for Gaussian and exponential spectrum are used to validate the proposed methods. For the case of unknown undulation height spectra of rough surfaces, the ACF is statistically calculated first, and then the KAC and KAW methods are used to calculate the average scattering coefficients numerically. Moreover, the numerical results show that the KAW method takes up more memory to retain enough energy of the ACF for accurate scattering coefficients. The two new methods are effective in the range of $\theta _{i} $ from 0° to 45° and have broad application prospects for real scenarios.