Empirical mode decomposition (EMD) is an adaptive method for nonlinear, non-stationary signal analysis. However, the upper and lower envelopes fitted by cubic spline interpolation (CSI) may often occur overshoots. In this paper, a new envelope fitting method based on the flattest constrained interpolation is proposed. The proposed method effectively integrates the difference between extremes into the cost function, and applies a chaos particle swarm optimization method to optimize the derivatives of the interpolation nodes. The proposed method was tested on three different types of data: ascertain signal, random signals and real electrocardiogram signals. The experimental results show that: (1) The proposed flattest envelope effectively solves the overshoots caused by CSI method and the artificial bends caused by piecewise parabola interpolation (PPI) method. (2) The index of orthogonality of the intrinsic mode functions (IMFs) based on the proposed method is 0.04054, 0.02222 ± 0.01468 and 0.04013 ± 0.03953 for the ascertain signal, random signals and electrocardiogram signals, respectively, which is lower than the CSI method and the PPI method, and means the IMFs are more orthogonal. (3) The index of energy conversation of the IMFs based on the proposed method is 0.96193, 0.93501 ± 0.03290 and 0.93041 ± 0.00429 for the ascertain signal, random signals and electrocardiogram signals, respectively, which is closer to 1 than the other two methods and indicates the total energy deviation amongst the components is smaller. (4) The comparisons of the Hilbert spectrums show that the proposed method overcomes the mode mixing problems very well, and make the instantaneous frequency more physically meaningful.