We present a rigorous analytical method for harmonic analysis of the angular error of rotary and linear encoders with sine/cosine output signals in quadrature that are distorted by superimposed Fourier series. To calculate the angle from measured sine and cosine encoder channels in quadrature, the arctangent function is commonly used. The hence non-linear relation between raw signals and calculated angle -- often thought of as a black box -- complicates the estimation of the angular error and its harmonic decomposition. By means of a Taylor series expansion of the harmonic amplitudes, our method allows for quantification of the impact of harmonic signal distortions on the angular error in terms of harmonic order, magnitude and phase, including an upper bound on the remaining error term -- without numerical evaluation of the arctangent function. The same approximation is achieved with an intuitive geometric approximation in the complex plane, validating the results. Additionally, interaction effects between harmonics in the signals are considered by higher-order Taylor expansion. The approximations show an excellent agreement with the exact calculation in numerical examples even in case of large distortion amplitudes, leading to practicable estimates for the angular error decomposition.