Python scipy.hanning() Examples

def ltsd_vad(x, fs, threshold=9, winsize=8192):
    # winsize based on sample rate
    # 1024 for fs = 16000
    orig_dtype = x.dtype
    orig_scale_min = x.min()
    orig_scale_max = x.max()
    x = (x - x.min()) / (x.max() - x.min())
    # works with 16 bit
    x = x * (2 ** 15)
    x = x.astype("int32")
    window = sp.hanning(winsize)
    ltsd = LTSD(winsize, window, 5)
    s_vad = ltsd.compute(x)
    # LTSD is 50% overlap, so each "step" covers 4096 samples
    # +1 to cover the extra edge window
    n_samples = int(((len(s_vad) + 1) * winsize) // 2)
    time_s = n_samples / float(fs)
    time_points = np.linspace(0, time_s, len(s_vad))
    time_samples = (fs * time_points).astype(np.int32)
    time_samples = time_samples
    f_vad = np.zeros_like(x, dtype=np.bool)
    offset = winsize
    for n, (ss, es) in enumerate(zip(time_samples[:-1], time_samples[1:])):
        sss = ss - offset
        if sss < 0:
            sss = 0
        ses = es - offset
        if ses < 0:
            ses = 0
        if s_vad[n + 1] < threshold:
            f_vad[sss:ses] = False
        else:
            f_vad[sss:ses] = True
    f_vad[ses:] = False
    x = x.astype("float64")
    x = x / float(2 ** 15)
    x = x * (orig_scale_max - orig_scale_min) + orig_scale_min
    x = x.astype(orig_dtype)
    return x[f_vad], f_vad 

def stft(x, fs, framesz, hop):
    framesamp = int(framesz*fs)
    hopsamp = int(hop*fs)
    w = scipy.hanning(framesamp)
    X = scipy.array([scipy.fft(w*x[i:i+framesamp])
                     for i in range(0, len(x)-framesamp, hopsamp)])
    return X 

def create_window(size):
    """
    A normalized Hanning window of given size. Useful for analyzing sinusoidal
    signals.
    """
    return normalized_window(scipy.hanning(size)) 

def test_hanning():
    """ Compare scipy and Matlab hanning window.
        Matlab returns a N+2 size window without first and last samples"""
    hanning = scipy.hanning(N_FRAME+2)[1:-1]
    hanning_m = eng.hanning(float(N_FRAME))
    hanning_m = np.array(hanning_m._data)
    assert_allclose(hanning, hanning_m, atol=ATOL) 

def test_hanning():
    """ Compare scipy and Matlab hanning window.

        Matlab returns a N+2 size window without first and last samples.
        A custom Octave function has been written to mimic this
        behavior."""
    hanning = scipy.hanning(N_FRAME+2)[1:-1]
    hanning_m = np.squeeze(octave.feval('octave/ml_hanning.m', N_FRAME))
    assert_allclose(hanning, hanning_m, atol=ATOL) 

def ltsd_vad(x, fs, threshold=9, winsize=8192):
    # winsize based on sample rate
    # 1024 for fs = 16000
    orig_dtype = x.dtype
    orig_scale_min = x.min()
    orig_scale_max = x.max()
    x = (x - x.min()) / (x.max() - x.min())
    # works with 16 bit
    x = x * (2 ** 15)
    x = x.astype("int32")
    window = sp.hanning(winsize)
    ltsd = LTSD(winsize, window, 5)
    s_vad = ltsd.compute(x)
    # LTSD is 50% overlap, so each "step" covers 4096 samples
    # +1 to cover the extra edge window
    n_samples = int(((len(s_vad) + 1) * winsize) // 2)
    time_s = n_samples / float(fs)
    time_points = np.linspace(0, time_s, len(s_vad))
    time_samples = (fs * time_points).astype(np.int32)
    time_samples = time_samples
    f_vad = np.zeros_like(x, dtype=np.bool)
    offset = winsize
    for n, (ss, es) in enumerate(zip(time_samples[:-1], time_samples[1:])):
        sss = ss - offset
        if sss < 0:
            sss = 0
        ses = es - offset
        if ses < 0:
            ses = 0
        if s_vad[n + 1] < threshold:
            f_vad[sss:ses] = False
        else:
            f_vad[sss:ses] = True
    f_vad[ses:] = False
    x = x.astype("float64")
    x = x / float(2 ** 15)
    x = x * (orig_scale_max - orig_scale_min) + orig_scale_min
    x = x.astype(orig_dtype)
    return x[f_vad], f_vad 

def ltsd_vad(x, fs, threshold=9, winsize=8192):
    # winsize based on sample rate
    # 1024 for fs = 16000
    orig_dtype = x.dtype
    orig_scale_min = x.min()
    orig_scale_max = x.max()
    x = (x - x.min()) / (x.max() - x.min())
    # works with 16 bit
    x = x * (2 ** 15)
    x = x.astype("int32")
    window = sp.hanning(winsize)
    ltsd = LTSD(winsize, window, 5)
    s_vad = ltsd.compute(x)
    # LTSD is 50% overlap, so each "step" covers 4096 samples
    # +1 to cover the extra edge window
    n_samples = int(((len(s_vad) + 1) * winsize) // 2)
    time_s = n_samples / float(fs)
    time_points = np.linspace(0, time_s, len(s_vad))
    time_samples = (fs * time_points).astype(np.int32)
    time_samples = time_samples
    f_vad = np.zeros_like(x, dtype=np.bool)
    offset = winsize
    for n, (ss, es) in enumerate(zip(time_samples[:-1], time_samples[1:])):
        sss = ss - offset
        if sss < 0:
            sss = 0
        ses = es - offset
        if ses < 0:
            ses = 0
        if s_vad[n + 1] < threshold:
            f_vad[sss:ses] = False
        else:
            f_vad[sss:ses] = True
    f_vad[ses:] = False
    x = x.astype("float64")
    x = x / float(2 ** 15)
    x = x * (orig_scale_max - orig_scale_min) + orig_scale_min
    x = x.astype(orig_dtype)
    return x[f_vad], f_vad 

def sinusoid_analysis(X, input_sample_rate, resample_block=128, copy=True):
    """
    Contruct a sinusoidal model for the input signal.

    Parameters
    ----------
    X : ndarray
        Input signal to model

    input_sample_rate : int
        The sample rate of the input signal

    resample_block : int, optional (default=128)
       Controls the step size of the sinusoidal model

    Returns
    -------
    frequencies_hz : ndarray
       Frequencies for the sinusoids, in Hz.

    magnitudes : ndarray
       Magnitudes of sinusoids returned in ``frequencies``

    References
    ----------
    D. P. W. Ellis (2004), "Sinewave Speech Analysis/Synthesis in Matlab",
    Web resource, available: http://www.ee.columbia.edu/ln/labrosa/matlab/sws/
    """
    X = np.array(X, copy=copy)
    resample_to = 8000
    if input_sample_rate != resample_to:
        if input_sample_rate % resample_to != 0:
            raise ValueError("Input sample rate must be a multiple of 8k!")
        # Should be able to use resample... ?
        # resampled_count = round(len(X) * resample_to / input_sample_rate)
        # X = sg.resample(X, resampled_count, window=sg.hanning(len(X)))
        X = sg.decimate(X, input_sample_rate // resample_to, zero_phase=True)
    step_size = 2 * round(resample_block / input_sample_rate * resample_to / 2.)
    a, g, e = lpc_analysis(X, order=8, window_step=step_size,
                           window_size=2 * step_size)
    f, m = lpc_to_frequency(a, g)
    f_hz = f * resample_to / (2 * np.pi)
    return f_hz, m 

def __init__(self, signal_frames, window=scipy.hanning, positive_only=True):
        """
        :param signal_frames: signal represented as SignalFrames instance
        :param window: STFT window function - produces 1D window which will
        be normalized
        """
        self.signal_frames = signal_frames

        x_frames = signal_frames.frames
        w = normalized_window(window(signal_frames.frame_size))

        # complex spectra of windowed blocks of signal - STFT
        self.X_complex = np.fft.fft(x_frames * w)
        # linear magnitude spectrogram
        self.X_mag = abs(self.X_complex) / self.X_complex.shape[1]

        # spectra of signal shifted in time

        # This fakes looking at the previous frame shifted by one sample.
        # In order to work only with one frame of size N and not N + 1, we fill the
        # missing value with zero. This should not introduce a large error, since the
        # borders of the amplitude frame will go to zero anyway due to applying a
        # window function in the STFT tranform.
        X_prev_time = np.fft.fft(shift_right(x_frames) * w)

        # spectra shifted in frequency
        X_prev_freq = shift_right(self.X_complex)

        # cross-spectra - ie. spectra of cross-correlation between the
        # respective time-domain signals
        X_cross_time = cross_spectrum(self.X_complex, X_prev_time)
        X_cross_freq = cross_spectrum(self.X_complex, X_prev_freq)

        # instantaneous frequency estimates
        # normalized frequencies in range [0.0, 1.0] - from DC to sample rate
        self.X_inst_freqs = estimate_instant_freqs(X_cross_time)
        # instantaneous group delay estimates
        # relative coordinates within the frame with range [-0.5, 0.5] where
        # 0.0 is the frame center
        self.X_group_delays = estimate_group_delays(X_cross_freq)

        if positive_only:
            self.X_mag = positive_freq_magnitudes(self.X_mag)
            self.X_complex, self.X_inst_freqs, self.X_group_delays = [
                select_positive_freq_fft(values) for values in
                [self.X_complex, self.X_inst_freqs, self.X_group_delays]
            ] 

def sinusoid_analysis(X, input_sample_rate, resample_block=128, copy=True):
    """
    Contruct a sinusoidal model for the input signal.

    Parameters
    ----------
    X : ndarray
        Input signal to model

    input_sample_rate : int
        The sample rate of the input signal

    resample_block : int, optional (default=128)
       Controls the step size of the sinusoidal model

    Returns
    -------
    frequencies_hz : ndarray
       Frequencies for the sinusoids, in Hz.

    magnitudes : ndarray
       Magnitudes of sinusoids returned in ``frequencies``

    References
    ----------
    D. P. W. Ellis (2004), "Sinewave Speech Analysis/Synthesis in Matlab",
    Web resource, available: http://www.ee.columbia.edu/ln/labrosa/matlab/sws/
    """
    X = np.array(X, copy=copy)
    resample_to = 8000
    if input_sample_rate != resample_to:
        if input_sample_rate % resample_to != 0:
            raise ValueError("Input sample rate must be a multiple of 8k!")
        # Should be able to use resample... ?
        # resampled_count = round(len(X) * resample_to / input_sample_rate)
        # X = sg.resample(X, resampled_count, window=sg.hanning(len(X)))
        X = sg.decimate(X, input_sample_rate // resample_to, zero_phase=True)
    step_size = 2 * round(resample_block / input_sample_rate * resample_to / 2.)
    a, g, e = lpc_analysis(X, order=8, window_step=step_size,
                           window_size=2 * step_size)
    f, m = lpc_to_frequency(a, g)
    f_hz = f * resample_to / (2 * np.pi)
    return f_hz, m 

def sinusoid_analysis(X, input_sample_rate, resample_block=128, copy=True):
    """
    Contruct a sinusoidal model for the input signal.

    Parameters
    ----------
    X : ndarray
        Input signal to model

    input_sample_rate : int
        The sample rate of the input signal

    resample_block : int, optional (default=128)
       Controls the step size of the sinusoidal model

    Returns
    -------
    frequencies_hz : ndarray
       Frequencies for the sinusoids, in Hz.

    magnitudes : ndarray
       Magnitudes of sinusoids returned in ``frequencies``

    References
    ----------
    D. P. W. Ellis (2004), "Sinewave Speech Analysis/Synthesis in Matlab",
    Web resource, available: http://www.ee.columbia.edu/ln/labrosa/matlab/sws/
    """
    X = np.array(X, copy=copy)
    resample_to = 8000
    if input_sample_rate != resample_to:
        if input_sample_rate % resample_to != 0:
            raise ValueError("Input sample rate must be a multiple of 8k!")
        # Should be able to use resample... ?
        # resampled_count = round(len(X) * resample_to / input_sample_rate)
        # X = sg.resample(X, resampled_count, window=sg.hanning(len(X)))
        X = sg.decimate(X, input_sample_rate // resample_to, zero_phase=True)
    step_size = 2 * round(resample_block / input_sample_rate * resample_to / 2.)
    a, g, e = lpc_analysis(X, order=8, window_step=step_size,
                           window_size=2 * step_size)
    f, m = lpc_to_frequency(a, g)
    f_hz = f * resample_to / (2 * np.pi)
    return f_hz, m