Ensemble empirical mode decomposition (EEMD) represents a valuable aid in empirical mode decomposition (EMD) and has been widely used in fault diagnosis of rolling element bearings. However, the intrinsic mode functions (IMFs) generated by EEMD often contain residual noise. In addition, adding different white Gaussian noise to the signal to be analyzed probably produces a different number of IMFs, and different number of IMFs makes difficult the averaging. To alleviate these two drawbacks, complete ensemble empirical mode decomposition with adaptive noise (CEEMDAN) was previously presented. Utilizing the advantages of CEEMDAN in extracting weak characteristics from noisy signals, a new fault diagnosis method of rolling element bearings based on CEEMDAN is proposed. With this method, a particular noise is added at each stage and after each IMF extraction, a unique residue is computed. In this way, this method solves the problem of the final averaging and obtains IMFs with less noise. A simulated signal is used to illustrate the effectiveness of the proposed method, and the decomposition results show that the method obtains more accurate IMFs than the EEMD. To further demonstrate the proposed method, it is applied to fault diagnosis of locomotive rolling element bearings. The diagnosis results prove that the method based on CEEMDAN may reveal the fault characteristic information of rolling element bearings better.