我们从Python开源项目中,提取了以下21个代码示例,用于说明如何使用scipy.signal.hamming()。
def __init__(self, win_len=0.025, win_step=0.01, num_filt=40, nfft=512, low_freq=20, high_freq=7800, pre_emph=0.97, win_fun=signal.hamming, **kwargs): super(FBank, self).__init__(**kwargs) if high_freq > self.fs / 2: raise ValueError("high_freq must be less or equal than fs/2") self.win_len = win_len self.win_step = win_step self.num_filt = num_filt self.nfft = nfft self.low_freq = low_freq self.high_freq = high_freq or self.fs / 2 self.pre_emph = pre_emph self.win_fun = win_fun self._filterbanks = self._get_filterbanks() self._num_feats = self.num_filt
def single_spectrogram(inseq,fs,wlen,h,imag=False): """ imag: Return Imaginary Data of the STFT on True """ NFFT = int(2**(np.ceil(np.log2(wlen)))) K = np.sum(hamming(wlen, False))/wlen raw_data = inseq.astype('float32') raw_data = raw_data/np.amax(np.absolute(raw_data)) stft_data,_,_ = STFT(raw_data,wlen,h,NFFT,fs) s = np.absolute(stft_data)/wlen/K; if np.fmod(NFFT,2): s[1:,:] *=2 else: s[1:-2] *=2 real_data = np.transpose(20*np.log10(s + 10**-6)).astype(np.float32) if imag: imag_data = np.angle(stft_data).astype(np.float32) return real_data,imag_data return real_data
def GLA(wsz, hop, N=4096): """ LSEE-MSTFT algorithm for computing the synthesis window used in inverse STFT method below. Args: wsz : (int) Synthesis window size hop : (int) Hop size N : (int) DFT Size Returns : symw: (array) Synthesised windowing function References : [1] Daniel W. Griffin and Jae S. Lim, ``Signal estimation from modified short-time Fourier transform,'' IEEE Transactions on Acoustics, Speech and Signal Processing, vol. 32, no. 2, pp. 236-243, Apr 1984. """ synw = hamming(wsz)/np.sqrt(N) synwProd = synw ** 2. synwProd.shape = (wsz, 1) redundancy = wsz/hop env = np.zeros((wsz, 1)) for k in xrange(-redundancy, redundancy + 1): envInd = (hop*k) winInd = np.arange(1, wsz+1) envInd += winInd valid = np.where((envInd > 0) & (envInd <= wsz)) envInd = envInd[valid] - 1 winInd = winInd[valid] - 1 env[envInd] += synwProd[winInd] synw = synw/env[:, 0] return synw
def amplitude_stream(data, sr, fftn, step, lowcut, highcut): '''returns an iterator with the start time in seconds and threshold value for each chunk''' from scipy.signal import hamming all_fft_freqs = np.fft.fftfreq(fftn, sr** -1) fft_freqs = (all_fft_freqs >= lowcut) & (all_fft_freqs <= highcut) window = hamming(fftn) data = data.ravel() for i in range(0, len(data) - fftn, step): x = data[i:i + fftn] * window fft = np.fft.fft(x, fftn)[fft_freqs] time = (i + fftn / 2) / sr yield time, np.mean(np.log(np.abs(fft)))
def __init__(self): self.outputNames=("OUTPUT",) # one output terminal named "output" self.outputTypes=("prim_int",) # the type is primitive int. self.data_array = [] self.prefix = './temp' self.wlen = 512 self.K = np.sum(hamming(self.wlen, False))/self.wlen
def wav_to_image(filename, wlen, mindata, maxdata, save=False, name_save=None, ): h = wlen/4 K = np.sum(hamming(wlen, False))/wlen nfft = int(2**(np.ceil(np.log2(wlen)))) Fs, data_seq = wavfile.read(filename) raw_data = data_seq.astype('float32') max_dt = np.amax(np.absolute(raw_data)) raw_data = raw_data/max_dt stft_data,_,_ = STFT(raw_data,wlen,h,nfft,Fs) s = abs(stft_data)/wlen/K; if np.fmod(nfft,2): s[1:,:] *=2 else: s[1:-2] *=2 data_temp = 20*np.log10(s + 10**-6) outdata = data_temp.transpose() """Scaling""" mindata = np.amin(outdata, axis=0, keepdims = True) maxdata = np.amax(outdata, axis=0, keepdims = True) outdata -=mindata outdata /=(maxdata-mindata) outdata *=0.8 outdata +=0.1 figmin = np.zeros((5,outdata.shape[1])) figmax = np.ones((5,outdata.shape[1])) outdata = np.concatenate((outdata,figmin,figmax), axis=0) dpi = 96 a = float(outdata.shape[0])/dpi b = float(outdata.shape[1])/dpi f = plt.figure(figsize=(b,a), dpi=dpi) f.figimage(outdata) if save: f.savefig(name_save, dpi=f.dpi) return f
def ISTFT(data, h, nfft, fs): # function: [x, t] = istft(stft, h, nfft, fs) # stft - STFT matrix (only unique points, time across columns, freq across rows) # h - hop size # nfft - number of FFT points # fs - sampling frequency, Hz # x - signal in the time domain # t - time vector, s # estimate the length of the signal coln = data.shape[1] xlen = nfft + (coln-1)*h x = np.zeros((xlen,)) # form a periodic hamming window win = hamming(nfft, False) # perform IFFT and weighted-OLA if np.fmod(nfft,2): lst_idx = -1 else: lst_idx = -2 for b in range (0, h*(coln-1),h): # extract FFT points X = data[:,1+b/h] X = np.concatenate((X, np.conjugate(X[lst_idx:0:-1]))) # IFFT xprim = np.real(np.fft.ifft(X)) # weighted-OLA x[b:b+nfft] = x[b:b+nfft] + np.transpose(xprim*win) W0 = np.sum(win*win) x *= h/W0 # calculate the time vector actxlen = x.shape[0] t = np.arange(0,actxlen-1,dtype=np.float32)/fs return x, t
def GLA(wsz, hop): """ LSEE-MSTFT algorithm for computing the synthesis window used in inverse STFT method below. Args: wsz : (int) Synthesis Window size hop : (int) Hop size Returns : symw: (array) Synthesised time-domain real signal. References : [1] Daniel W. Griffin and Jae S. Lim, ``Signal estimation from modified short-time Fourier transform,'' IEEE Transactions on Acoustics, Speech and Signal Processing, vol. 32, no. 2, pp. 236-243, Apr 1984. """ synw = hamming(wsz)/np.sum(hamming(wsz)) synwProd = synw ** 2. synwProd.shape = (wsz, 1) redundancy = wsz/hop env = np.zeros((wsz, 1)) for k in xrange(-redundancy, redundancy + 1): envInd = (hop*k) winInd = np.arange(1, wsz+1) envInd += winInd valid = np.where((envInd > 0) & (envInd <= wsz)) envInd = envInd[valid] - 1 winInd = winInd[valid] - 1 env[envInd] += synwProd[winInd] synw = synw/env[:, 0] return synw
def get_ft_windows(): """ Retrieve the available windows to be applied on signal data before FT. @return: dict with keys being the window name and items being again a dict containing the actual function and the normalization factor to calculate correctly the amplitude spectrum in the Fourier Transform To find out the amplitude normalization factor check either the scipy implementation on https://github.com/scipy/scipy/blob/v0.15.1/scipy/signal/windows.py#L336 or just perform a sum of the window (oscillating parts of the window should be averaged out and constant offset factor will remain): MM=1000000 # choose a big number print(sum(signal.hanning(MM))/MM) """ win = {'none': {'func': np.ones, 'ampl_norm': 1.0}, 'hamming': {'func': signal.hamming, 'ampl_norm': 1.0/0.54}, 'hann': {'func': signal.hann, 'ampl_norm': 1.0/0.5}, 'blackman': {'func': signal.blackman, 'ampl_norm': 1.0/0.42}, 'triang': {'func': signal.triang, 'ampl_norm': 1.0/0.5}, 'flattop': {'func': signal.flattop, 'ampl_norm': 1.0/0.2156}, 'bartlett': {'func': signal.bartlett, 'ampl_norm': 1.0/0.5}, 'parzen': {'func': signal.parzen, 'ampl_norm': 1.0/0.375}, 'bohman': {'func': signal.bohman, 'ampl_norm': 1.0/0.4052847}, 'blackmanharris': {'func': signal.blackmanharris, 'ampl_norm': 1.0/0.35875}, 'nuttall': {'func': signal.nuttall, 'ampl_norm': 1.0/0.3635819}, 'barthann': {'func': signal.barthann, 'ampl_norm': 1.0/0.5}} return win
def convert_input(channel, annotation): """ Into output """ # Remove non-beat annotations beats = beat_annotations(annotation) # Create dirac-comb signal dirac = np.zeros_like(channel) dirac[beats] = 1.0 # Use hamming window as a bell-curve filter width = 36 filter = ss.hamming(width) gauss = np.convolve(filter, dirac, mode = 'same') return dirac, gauss
def CQT_fast(x,fs,bins,fmin,fmax,M): threshold = 0.0054 #for Hamming window K = int(bins*np.ceil(np.log2(fmax/fmin))) Q = 1/(2**(1/bins)-1) nfft = np.int32(nearestPow2(np.ceil(Q*fs/fmin))) tempKernel = np.zeros(nfft, dtype = np.complex) specKernel = np.zeros(nfft, dtype = np.complex) sparKernel = [] #create sparse Kernel for k in range(K-1,-1,-1): fk = (2**(k/bins))*fmin N = np.int32(np.round((Q*fs)/fk)) tempKernel[:N] = hamming(N)/N * np.exp(-2*np.pi*1j*Q*np.arange(N)/N) specKernel = fft(tempKernel) specKernel[np.where(np.abs(specKernel) <= threshold)] = 0 if k == K-1: sparKernel = specKernel else: sparKernel = np.vstack((specKernel, sparKernel)) sparKernel = np.transpose(np.conjugate(sparKernel))/nfft ft = fft(x,nfft) cqt = np.dot(ft, sparKernel) ft = fft(x,nfft*(2**M)) #calculate harmonic power spectrum #harm_pow = HPS(ft,M) #cqt = np.dot(harm_pow, sparKernel) return cqt
def CQT_slow(x, fs, bins, fmin, fmax): K = int(bins*np.ceil(np.log2(fmax/fmin))) Q = 1/(2**(1/bins)-1) cqt = np.zeros(K, dtype = np.complex) for k in range(K): fk = (2**(k/bins))*fmin N = int(np.round(Q*fs/fk)) arr = -2*np.pi*1j*Q*np.arange(N)/N cqt[k] = np.dot(x[:N], np.transpose(hamming(N) * np.exp(arr)))/N return cqt
def windowing(input): """ Applies hamming window to the input frames. Args: input: array of speech samples [N x M] where N is the number of frames and M the samples per frame Output: array of windoed speech samples [N x M] Note (you can use the function hamming from scipy.signal, include the sym=0 option if you want to get the same results as in the example) """ N, M = input.shape window = ss.hamming(M, sym=False) return (input * window)
def mfcc(s,fs, nfiltbank): #divide into segments of 25 ms with overlap of 10ms nSamples = np.int32(0.025*fs) overlap = np.int32(0.01*fs) nFrames = np.int32(np.ceil(len(s)/(nSamples-overlap))) #zero padding to make signal length long enough to have nFrames padding = ((nSamples-overlap)*nFrames) - len(s) if padding > 0: signal = np.append(s, np.zeros(padding)) else: signal = s segment = np.empty((nSamples, nFrames)) start = 0 for i in range(nFrames): segment[:,i] = signal[start:start+nSamples] start = (nSamples-overlap)*i #compute periodogram nfft = 512 periodogram = np.empty((nFrames,nfft/2 + 1)) for i in range(nFrames): x = segment[:,i] * hamming(nSamples) spectrum = fftshift(fft(x,nfft)) periodogram[i,:] = abs(spectrum[nfft/2-1:])/nSamples #calculating mfccs fbank = mel_filterbank(nfft, nfiltbank, fs) #nfiltbank MFCCs for each frame mel_coeff = np.empty((nfiltbank,nFrames)) for i in range(nfiltbank): for k in range(nFrames): mel_coeff[i,k] = np.sum(periodogram[k,:]*fbank[:,i]) mel_coeff = np.log10(mel_coeff) mel_coeff = dct(mel_coeff) #exclude 0th order coefficient (much larger than others) mel_coeff[0,:]= np.zeros(nFrames) return mel_coeff
def initLanczos(self, filtOrder): self.filtOrder = filtOrder if self.filtOrder % 2 != 0: raise Exception('Invalid filtOrder: ' + str(self.filtOrder) + ' Must be an even integer.') radius = self.filtOrder // 2 win = np.sinc(np.linspace(-radius, radius, self.filtOrder+1) / float(radius)) # lanczos #win = spsig.hamming(self.filtOrder+1) # sinc-hamming # this should be automated somehow XXX - idfah if self.filtOrder <= 6: cutoff = 2*0.570 elif self.filtOrder <= 8: cutoff = 2*0.676 elif self.filtOrder <= 12: cutoff = 2*0.781 elif self.filtOrder <= 16: cutoff = 2*0.836 elif self.filtOrder <= 32: cutoff = 2*0.918 elif self.filtOrder <= 64: cutoff = 2*0.959 # need to fix for multiple pool sizes XXX - idfah cutoff /= float(self.poolSize) taps = cutoff * np.linspace(-radius, radius, self.filtOrder+1, dtype=self.dtype) impulseResponse = cutoff * np.sinc(taps) * win self.filters = [] nReadoutLayers = 1 if self.nHidden is None else 2 for ni, no in self.layerDims[:-nReadoutLayers]: noEmb = no*(self.filtOrder+1) # no outs after filter embedding filtMat = np.zeros(noEmb*2, dtype=self.dtype) filtMat[noEmb-1::-no] = impulseResponse # filters strided for embedding sz = filtMat.itemsize filtMat = npst.as_strided(filtMat, (no,noEmb), strides=(sz,sz))[::-1].T self.filters.append(filtMat.copy())
def STFT(x, wlen, h, nfft, fs): ######################################################## # Short-Time Fourier Transform % # with MATLAB Implementation % # For Python % # Copier: Nelson Yalta 11/03/15 % ######################################################## # function: [stft, f, t] = stft(x, wlen, h, nfft, fs) # x - signal in the time domain # wlen - length of the hamming window # h - hop size # nfft - number of FFT points # fs - sampling frequency, Hz # f - frequency vector, Hz # t - time vector, s # stft - STFT matrix (only unique points, time across columns, freq across rows) # represent x as column-vector if it is not if (len(x.shape) > 1) and (x.shape[1] > 1): x = x.transpose() # length of the signal xlen = x.shape[0] # form a periodic hamming window win = hamming(wlen, False) # form the stft matrix rown = int(np.ceil((1.0+nfft)/2)) coln = int(np.fix((xlen-wlen)/h) + 1) short_tft = np.zeros((rown,coln)).astype('complex64') # initialize the indexes indx = 0 col = 0 # perform STFT while (indx + wlen <= xlen): # windowing xw =x[indx:indx+wlen]*win # FFT X = np.fft.fft(xw,nfft) # update the stft matrix short_tft[:,col] = X[0:rown] # update the indexes indx += h col += 1 # calculate the time and frequency vectors t = np.linspace(wlen/2,wlen/2+(coln-1)*h,coln)/fs f = np.arange(0,rown,dtype= np.float32)*fs/nfft return short_tft, f, t
def show_objective(): """ For the model """ # Choose a record records = dm.get_records() path = records[17] record = wf.rdsamp(path) ann = wf.rdann(path, 'atr') chid = 0 print 'Channel:', record.signame[chid] cha = record.p_signals[:, chid] # These were found manually sta = 184000 end = sta + 1000 times = np.arange(end-sta, dtype = 'float') times /= record.fs # Extract the annotations for that fragment where = (sta < ann.annsamp) & (ann.annsamp < end) samples = ann.annsamp[where] - sta print samples # Prepare dirac-comb type of labels qrs_values = np.zeros_like(times) qrs_values[samples] = 1 # Prepare gaussian-comb type of labels kernel = ss.hamming(36) qrs_gauss = np.convolve(kernel, qrs_values, mode = 'same') # Make the plots fig = plt.figure() ax1 = fig.add_subplot(3,1,1) ax1.plot(times, cha[sta : end]) ax2 = fig.add_subplot(3,1,2, sharex=ax1) ax2.plot(times, qrs_values, 'C1', lw = 4, alpha = 0.888) ax3 = fig.add_subplot(3,1,3, sharex=ax1) ax3.plot(times, qrs_gauss, 'C3', lw = 4, alpha = 0.888) plt.setp(ax1.get_xticklabels(), visible=False) plt.setp(ax2.get_xticklabels(), visible=False) plt.xlabel('Time [s]') plt.xlim([0, 2.5]) plt.show()
def show_objective_part2(): """ For the model """ # Choose a record records = dm.get_records() path = records[13] record = wf.rdsamp(path) ann = wf.rdann(path, 'atr') chid = 0 print 'File:', path print 'Channel:', record.signame[chid] cha = record.p_signals[:, chid] # These were found manually sta = 184000 end = sta + 1000 times = np.arange(end-sta, dtype = 'float') times /= record.fs # Extract the annotations for that fragment where = (sta < ann.annsamp) & (ann.annsamp < end) samples = ann.annsamp[where] - sta print samples # Prepare dirac-comb type of labels qrs_values = np.zeros_like(times) qrs_values[samples] = 1 # Prepare gaussian-comb type of labels kernel = ss.hamming(36) qrs_gauss = np.convolve(kernel, qrs_values, mode = 'same') # Make the plots fig = plt.figure() ax1 = fig.add_subplot(2,1,1) ax1.plot(times, cha[sta : end]) ax1.set_title('Input', loc = 'left') ax2 = fig.add_subplot(2,1,2, sharex=ax1) ax2.plot(times, qrs_gauss, 'C3', lw = 4, alpha = 0.888) ax2.set_title('Output', loc = 'left') ax1.grid() ax2.grid() plt.setp(ax1.get_xticklabels(), visible=False) plt.xlabel('Time [s]') plt.xlim([0, 2.5]) plt.show()