In [2]:

    

from matplotlib.mlab import psd

def hammingwindow(freq1, freq2, f_sample=10, width=3, sig=.5):
    "Erstellt einen hamming windowed Frequenzübergang Freq1 auf Freq2 mit N Perioden Freq1 Rechteckfunktion+selbe Anz Data Freq2 und entsprechende FFT"
    freq1, freq2, width = float(freq1), float(freq2), float(width)
    f_sample = f_sample*1000
    #print "Samplingfrequenz:", f_sample/1000, "kHz"
    f1 = tile(concatenate((ones(f_sample/(2*freq1)), zeros(f_sample/(2*freq1))), axis=1), sig*freq1)
    f2 = tile(concatenate((ones(f_sample/(2*freq2)), zeros(f_sample/(2*freq2))), axis=1), sig*freq2)
    f = concatenate((f1, f2), axis=1)-0.5
    xaxis = linspace(0, 2*sig, f.size)
    window = concatenate((zeros(((width-1)/(2*width))*f.size), numpy.hamming(f.size/width)), axis=1)
    window = concatenate((window, zeros(((width-1)/(2*width))*f.size)), axis=1)
    window = concatenate((window, zeros(f.size-window.size)), axis=1)
    f = f*window
    #plot(xaxis, f)
    #plt.xlabel('t in s')
    #plt.ylim(-1,1)
    #figure()
    #FT = fft.fft(f)/f.size
    #FT = FT[:(FT.size/2)]
    #faxis = fft.fftfreq(f.size, (1./f_sample))
    #faxis = faxis[:(faxis.size/2)]
    #plot(faxis, abs(FT))
    #plt.xlim(0, (freq1+freq2))
    #plt.xlabel('Freq in Hz')
    #plt.ylabel('Amplitude')
    #figure()
    s, freqs = psd(f, 29929, f_sample)
    #plt.xlim(0, (freq1+freq2))
    return 20*log10(s), freqs;

def hammingwindow3d(freq1, freq2, f_sample=10, data=100, step=0.1):
    FT, faxis = hammingwindow(freq1, freq2, f_sample)
    print "Samplingfrequenz:", f_sample, "kHz"
    S = zeros((data-10, FT.size))
    l = arange(1, data*step, step)
    for n, item in enumerate(l):
        sn, faxis = hammingwindow(freq1, freq2, f_sample, item)
        S[n,:] = sn.T
    return S, faxis, l;

S, faxis, l = hammingwindow3d(300,800,48)
'''
from mpl_toolkits.mplot3d import Axes3D

x, y = meshgrid(faxis[:100], l)
fig = plt.figure()
ax = fig.add_subplot(1,1,1, projection='3d')
ax.plot_surface(x, y, S[:,:100], cmap=plt.get_cmap("autumn"))
plt.show()
'''

figure()
im = imshow(S[:,:1000], extent=[0, faxis[1000], 9, 0], aspect=('auto'))
plt.colorbar()

Smax = amax(S, axis=0)
S_norm = zeros(S.shape)
for i, item in enumerate(Smax):
    S_norm[:,i] = -item/abs(S[:,i])


figure()
im = imshow(S_norm[:,:1000], extent=[0, faxis[1000], 9, 0], aspect=('auto'), cmap='hot')
plt.colorbar()


    

    




Samplingfrequenz: 48 kHz






    
Out[2]:





<matplotlib.colorbar.Colorbar instance at 0x7f7c4e35b2d8>






    













    

In [ ]: