MATLAB实现变分模态分解(VMD)代码

function [u, u_hat, omega] = VMD(signal, alpha, tau, K, DC, init, tol)
% Variational Mode Decomposition
% Authors: Konstantin Dragomiretskiy and Dominique Zosso
% zosso@math.ucla.edu --- http://www.math.ucla.edu/~zosso
% Initial release 2013-12-12 (c) 2013
%
% Input and Parameters:
% ---------------------
% signal  - the time domain signal (1D) to be decomposed         
% alpha   - the balancing parameter of the data-fidelity         
%           constraint,一般设置为2000
% tau     - time-step of the dual ascent ( pick 0 for noise-slack )，       
% K       - the number of modes to be recovered，
% DC      - true if the first mode is put and kept at DC (0-freq)，
% init    - 0 = all omegas start at 0                               
%                    1 = all omegas start uniformly distributed，
%                    2 = all omegas initialized randomly，
% tol     - tolerance of convergence criterion; typically around 1e-6  
%
% Output:
% -------
% u       - the collection of decomposed modes              （分解的模态分量）
% u_hat   - spectra of the modes，
% omega   - estimated mode center-frequencies，
% When using this code, please do cite our paper:
% -----------------------------------------------
% K. Dragomiretskiy, D. Zosso, Variational Mode Decomposition, IEEE Trans.
% on Signal Processing (in press)
% please check here for update reference: 
%          http://dx.doi.org/10.1109/TSP.2013.2288675

%---------- Preparations

% Period and sampling frequency of input signal
save_T = length(signal);
fs = 1/save_T;

% extend the signal by mirroring
T = save_T;
f_mirror(1:T/2) = signal(T/2:-1:1);
f_mirror(T/2+1:3*T/2) = signal;
f_mirror(3*T/2+1:2*T) = signal(T:-1:T/2+1);
f = f_mirror;

% Time Domain 0 to T (of mirrored signal)
T = length(f);
t = (1:T)/T;

% Spectral Domain discretization
freqs = t-0.5-1/T;

% Maximum number of iterations (if not converged yet, then it won't anyway)
N = 500;

% For future generalizations: individual alpha for each mode
Alpha = alpha*ones(1,K);

% Construct and center f_hat