laszukdawid / PyEMD

Python implementation of Empirical Mode Decompoisition (EMD) method
https://pyemd.readthedocs.io/
Apache License 2.0
867 stars 224 forks source link

Fix CEEMDAN #45

Closed nescirem closed 5 years ago

nescirem commented 5 years ago

My python version is 3.6.4. I have tested EMD, EEMD and CEEMDAN to decomposite a piece of music, and I found the plotted Original signal is wrong if I use CEEMDAN as you can see in Fig.3, compared with Figs.1-2. The same issue is opened in #42 . In #42 and #44 one can find he will receive a wrong result if he uses res = S-np.sum(C_IMFs, axis=0) to calculate the residual due to the wrong Original signal. music_emd.png Fig.1 EMD music_eemd.png Fig.2 EEMD music_ceemdan.png Fig.3 CEEMDAN


I debugged it and found out why the Original signal is changed: the function ceemdan changed it. It may need to be scaled inside the function but it shouldn't affect the outside one. That is why my copy action worked. But it dosesn't slove the problem because I just changed the usage of it.


I don't think input signal should changed outside the function ceemdan and fixed it in a better way instead of just change the usage of it. For the reason why it works, one can refer to https://stackoverflow.com/questions/22054698/python-modifying-list-inside-a-function The result of CEEMDAN after I fixed it: music_ceemdan_fixed.png Fig.4 Fixed CEEMDAN

Attention: One may need to run python setup.py install again to update the function.


The code when I use CEEMDAN to dcoposite the music (from:http://www.music.helsinki.fi/tmt/opetus/uusmedia/esim/index-e.html):

import wave 
import numpy as np
import pylab as plt
from PyEMD import CEEMDAN

max_imf = -1

# Import music
filename= 'a2002011001-e02-8kHz.wav'
f = wave.open(filename,'rb')
params = f.getparams()
nchannels, sampwidth, framerate, nframes = params[:4]
strData = f.readframes(nframes)
waveData = np.fromstring(strData,dtype=np.int16)
waveData = waveData*1.0/(max(abs(waveData))) #Normalize
waveData = np.reshape(waveData,[nframes,nchannels])
f.close()

time = np.arange(0,nframes)*(1.0 / framerate)

t = time[:4000]
tMin = t[1]
tMax = t[-1]
S = waveData[:4000,0]

# Prepare and run CEEMDAN
ceemdan = CEEMDAN()

C_IMFs = ceemdan(S, t, max_imf)
imfNo  = C_IMFs.shape[0]

# Plot results in a grid
c = np.floor(np.sqrt(imfNo+2))
r = np.ceil((imfNo+2)/c)

plt.ioff()
plt.figure(figsize=(10,10))
plt.subplot(r,c,1)
plt.plot(t, S, 'r')
plt.xlim((tMin, tMax))
plt.title("Original signal")
for num in range(imfNo):
    plt.subplot(r,c,num+2)
    plt.plot(t, C_IMFs[num],'g')
    plt.xlim((tMin, tMax))
    plt.title("Imf "+str(num+1))

plt.show()
laszukdawid commented 5 years ago

Thank you for contributing.

Could you please update this pull request to a single commit (squeeze them) instead of 4?

Also, the only change that is required is in CEEMDAN.py on line 172. Please don't update .gitignore at this time.

It's also a good practise to update test cases so that they reflect what you are trying to fix. Any chance you could also add a test case for what you're fixing? It should be with other CEEMDAN test, i.e. in https://github.com/laszukdawid/PyEMD/blob/master/PyEMD/tests/test_ceemdan.py .

laszukdawid commented 5 years ago

Improved PR in #46