Python scipy.fftpack.ifft() Examples

The following are 30 code examples of scipy.fftpack.ifft(). You can vote up the ones you like or vote down the ones you don't like, and go to the original project or source file by following the links above each example. You may also want to check out all available functions/classes of the module scipy.fftpack , or try the search function .
Example #1
Source Project: simnibs   Author: simnibs   File: grid.py    License: GNU General Public License v3.0 7 votes vote down vote up
def quadrature_cc_1D(N):
    """ Computes the Clenshaw Curtis nodes and weights """
    N = np.int(N)        
    if N == 1:
        knots = 0
        weights = 2
    else:
        n = N - 1
        C = np.zeros((N,2))
        k = 2*(1+np.arange(np.floor(n/2)))
        C[::2,0] = 2/np.hstack((1, 1-k*k))
        C[1,1] = -n
        V = np.vstack((C,np.flipud(C[1:n,:])))
        F = np.real(ifft(V, n=None, axis=0))
        knots = F[0:N,1]
        weights = np.hstack((F[0,0],2*F[1:n,0],F[n,0]))
            
    return knots, weights 
Example #2
Source Project: ClearMap   Author: ChristophKirst   File: Convolution.py    License: GNU General Public License v3.0 6 votes vote down vote up
def _centered(arr, newsize):
    # Return the center newsize portion of the array.
    newsize = numpy.asarray(newsize)
    currsize = numpy.array(arr.shape)
    startind = (currsize - newsize) // 2
    endind = startind + newsize
    myslice = [slice(startind[k], endind[k]) for k in range(len(endind))]
    return arr[tuple(myslice)]


#def _rfftn(a, s=None, axes=None):
#    s, axes = _cook_nd_args(a, s, axes);
#    a = fft.rfft(a, s[-1], axes[-1], overwrite_x = True)
#    for ii in range(len(axes)-1):
#        a = fft.fft(a, s[ii], axes[ii], overwrite_x = True)
#    return a

#def _irfftn(a, s=None, axes=None):   
#    #a = asarray(a).astype('complex64')
#    s, axes = _cook_nd_args(a, s, axes, invreal=1)
#    for ii in range(len(axes)-1):
#        a = fft.ifft(a, s[ii], axes[ii], overwrite_x = True);
#    a = fft.ifft(a, s[-1], axes[-1], overwrite_x = True);
#    a = a.real;
#    return a 
Example #3
Source Project: dispel4py   Author: dispel4py   File: whiten.py    License: Apache License 2.0 6 votes vote down vote up
def spectralwhitening(stream):
    """
    Apply spectral whitening to data.
    Data is divided by its smoothed (Default: None) amplitude spectrum.
    """
    stream2 = copy.deepcopy(stream)

    for trace in arange(len(stream2)):
        data = stream2[trace].data

        n = len(data)
        nfft = nextpow2(n)

        spec = fft(data, nfft)
        spec_ampl = sqrt(abs(multiply(spec, conjugate(spec))))

        spec /= spec_ampl  # Do we need to do some smoothing here?
        ret = real(ifft(spec, nfft)[:n])

        stream2[trace].data = ret

    return stream2 
Example #4
Source Project: dispel4py   Author: dispel4py   File: whiten.py    License: Apache License 2.0 6 votes vote down vote up
def spectralwhitening_smooth(stream, N):
    """
    Apply spectral whitening to data.
    Data is divided by its smoothed (Default: None) amplitude spectrum.
    """
    stream2 = copy.deepcopy(stream)

    for trace in arange(len(stream2)):
        data = stream2[trace].data

        n = len(data)
        nfft = nextpow2(n)

        spec = fft(data, nfft)
        spec_ampl = sqrt(abs(multiply(spec, conjugate(spec))))

        spec_ampl = smooth(spec_ampl, N)

        spec /= spec_ampl  # Do we need to do some smoothing here?
        ret = real(ifft(spec, nfft)[:n])

        stream2[trace].data = ret

    return stream2 
Example #5
Source Project: fftoptionlib   Author: arraystream   File: fourier_pricer.py    License: BSD 3-Clause "New" or "Revised" License 6 votes vote down vote up
def carr_madan_fraction_fft_call_pricer(N, d_u, d_k, alpha, r, t, S0, q, chf_ln_st):
    rou = (d_u * d_k) / (2 * np.pi)
    beta = np.log(S0) - d_k * N / 2
    u_arr = np.arange(N) * d_u
    k_arr = beta + np.arange(N) * d_k
    delta_arr = np.zeros(N)
    delta_arr[0] = 1
    w_arr = d_u / 3 * (3 + (-1) ** (np.arange(N) + 1) - delta_arr)
    call_chf = (np.exp(-r * t) / ((alpha + 1j * u_arr) * (alpha + 1j * u_arr + 1))) * chf_ln_st(
        u_arr - (alpha + 1) * 1j,
        t, r, q=q, S0=S0)
    x_arr = np.exp(-1j * beta * u_arr) * call_chf * w_arr
    y_arr = np.zeros(2 * N) * 0j
    y_arr[:N] = np.exp(-1j * np.pi * rou * np.arange(N) ** 2) * x_arr
    z_arr = np.zeros(2 * N) * 0j
    z_arr[:N] = np.exp(1j * np.pi * rou * np.arange(N) ** 2)
    z_arr[N:] = np.exp(1j * np.pi * rou * np.arange(N - 1, -1, -1) ** 2)
    ffty = (fft(y_arr))
    fftz = (fft(z_arr))
    fftx = ffty * fftz
    fftpsi = ifft(fftx)
    fft_prices = np.exp(-1j * np.pi * (np.arange(N) ** 2) * rou) * fftpsi[:N]
    call_prices = (np.exp(-alpha * k_arr) / np.pi) * fft_prices.real
    return np.exp(k_arr), call_prices 
Example #6
Source Project: Computable   Author: ktraunmueller   File: test_basic.py    License: MIT License 6 votes vote down vote up
def test_definition(self):
        x1 = [1,2,3,4,1,2,3,4]
        x1_1 = [1,2+3j,4+1j,2+3j,4,2-3j,4-1j,2-3j]
        x2 = [1,2,3,4,1,2,3,4,5]
        x2_1 = [1,2+3j,4+1j,2+3j,4+5j,4-5j,2-3j,4-1j,2-3j]

        def _test(x, xr):
            y = irfft(np.array(x, dtype=self.rdt))
            y1 = direct_irdft(x)
            self.assertTrue(y.dtype == self.rdt,
                    "Output dtype is %s, expected %s" % (y.dtype, self.rdt))
            assert_array_almost_equal(y,y1, decimal=self.ndec)
            assert_array_almost_equal(y,ifft(xr), decimal=self.ndec)

        _test(x1, x1_1)
        _test(x2, x2_1) 
Example #7
Source Project: pwtools   Author: elcorto   File: test_pdos.py    License: BSD 3-Clause "New" or "Revised" License 6 votes vote down vote up
def test_pdos_1d():
    pad=lambda x: pad_zeros(x, nadd=len(x)-1)
    n=500; w=welch(n)
    # 1 second signal
    t=np.linspace(0,1,n); dt=t[1]-t[0]
    # sum of sin()s with random freq and phase shift, 10 frequencies from
    # f=0...100 Hz
    v=np.array([np.sin(2*np.pi*f*t + rand()*2*np.pi) for f in rand(10)*100]).sum(0)
    f=np.fft.fftfreq(2*n-1, dt)[:n]

    c1=mirror(ifft(abs(fft(pad(v)))**2.0)[:n].real)
    c2=correlate(v,v,'full')
    c3=mirror(acorr(v,norm=False))
    assert np.allclose(c1, c2)
    assert np.allclose(c1, c3)

    p1=(abs(fft(pad(v)))**2.0)[:n]
    p2=(abs(fft(mirror(acorr(v,norm=False)))))[:n]
    assert np.allclose(p1, p2)

    p1=(abs(fft(pad(v*w)))**2.0)[:n]
    p2=(abs(fft(mirror(acorr(v*w,norm=False)))))[:n]
    assert np.allclose(p1, p2) 
Example #8
Source Project: tools   Author: kastnerkyle   File: audio_tools.py    License: BSD 3-Clause "New" or "Revised" License 6 votes vote down vote up
def istft(X, fftsize=128, step="half", wsola=False, mean_normalize=True,
          real=False, compute_onesided=True):
    """
    Compute ISTFT for STFT transformed X
    """
    if real:
        local_ifft = fftpack.irfft
        X_pad = np.zeros((X.shape[0], X.shape[1] + 1)) + 0j
        X_pad[:, :-1] = X
        X = X_pad
    else:
        local_ifft = fftpack.ifft
    if compute_onesided:
        X_pad = np.zeros((X.shape[0], 2 * X.shape[1])) + 0j
        X_pad[:, :fftsize // 2 + 1] = X
        X_pad[:, fftsize // 2 + 1:] = 0
        X = X_pad
    X = local_ifft(X).astype("float64")
    if step == "half":
        X = invert_halfoverlap(X)
    else:
        X = overlap_add(X, step, wsola=wsola)
    if mean_normalize:
        X -= np.mean(X)
    return X 
Example #9
Source Project: GraphicDesignPatternByPython   Author: Relph1119   File: test_basic.py    License: MIT License 6 votes vote down vote up
def test_size_accuracy(self):
        # Sanity check for the accuracy for prime and non-prime sized inputs
        if self.rdt == np.float32:
            rtol = 1e-5
        elif self.rdt == np.float64:
            rtol = 1e-10

        for size in LARGE_COMPOSITE_SIZES + LARGE_PRIME_SIZES:
            np.random.seed(1234)
            x = np.random.rand(size).astype(self.rdt)
            y = ifft(fft(x))
            _assert_close_in_norm(x, y, rtol, size, self.rdt)
            y = fft(ifft(x))
            _assert_close_in_norm(x, y, rtol, size, self.rdt)

            x = (x + 1j*np.random.rand(size)).astype(self.cdt)
            y = ifft(fft(x))
            _assert_close_in_norm(x, y, rtol, size, self.rdt)
            y = fft(ifft(x))
            _assert_close_in_norm(x, y, rtol, size, self.rdt) 
Example #10
Source Project: fluids   Author: CalebBell   File: pychebfun.py    License: MIT License 6 votes vote down vote up
def polyval(self, chebcoeff):
        """
        Compute the interpolation values at Chebyshev points.
        chebcoeff: Chebyshev coefficients
        """
        N = len(chebcoeff)
        if N == 1:
            return chebcoeff

        data = even_data(chebcoeff)/2
        data[0] *= 2
        data[N-1] *= 2

        fftdata = 2*(N-1)*fftpack.ifft(data, axis=0)
        complex_values = fftdata[:N]
        # convert to real if input was real
        if np.isrealobj(chebcoeff):
            values = np.real(complex_values)
        else:
            values = complex_values
        return values 
Example #11
Source Project: pyGSTi   Author: pyGSTio   File: signal.py    License: Apache License 2.0 5 votes vote down vote up
def idft(modes, null_hypothesis, counts=1):
    """
    Inverts the dft function.

    Parameters
    ----------
    modes : array
        The fourier modes to be transformed to the time domain.

    null_hypothesis : array
        The array that was used in the normalization before the dct. This is
        commonly the mean of the time-domain data vector. All elements of this
        array must be in (0,1).

     counts : int, optional
        A factor in the normalization, that should correspond to the counts-per-timestep (so
        for full time resolution this is 1).

    Returns
    -------
    array
        Inverse of the dft function

    """
    z = _np.sqrt(len(modes)) * _ifft(modes)  # TIM CHECK THIS: len(*modes*) correct?
    x = unstandardizer(z, null_hypothesis, counts)
    return x 
Example #12
Source Project: nussl   Author: nussl   File: repet.py    License: MIT License 5 votes vote down vote up
def compute_beat_spectrum(power_spectrogram):
        """ Computes the beat spectrum averages (over freq's) the autocorrelation matrix of a one-sided spectrogram.

        The autocorrelation matrix is computed by taking the autocorrelation of each row of the spectrogram and
        dismissing the symmetric half.

        Args:
            power_spectrogram (:obj:`np.array`): 2D matrix containing the one-sided power spectrogram of an audio signal
            
        Returns:
            (:obj:`np.array`): array containing the beat spectrum based on the power spectrogram
            
        See Also:
            J Foote's original derivation of the Beat Spectrum: 
            Foote, Jonathan, and Shingo Uchihashi. "The beat spectrum: A new approach to rhythm analysis." 
            Multimedia and Expo, 2001. ICME 2001. IEEE International Conference on. IEEE, 2001.
            (`See PDF here <http://rotorbrain.com/foote/papers/icme2001.pdf>`_)
            
        """
        freq_bins, time_bins = power_spectrogram.shape

        # row-wise autocorrelation according to the Wiener-Khinchin theorem
        power_spectrogram = np.vstack([power_spectrogram, np.zeros_like(power_spectrogram)])

        nearest_power_of_two = 2 ** np.ceil(np.log(power_spectrogram.shape[0]) / np.log(2))
        pad_amount = int(nearest_power_of_two - power_spectrogram.shape[0])
        power_spectrogram = np.pad(power_spectrogram, ((0, pad_amount), (0, 0)), 'constant')
        fft_power_spec = scifft.fft(power_spectrogram, axis=0)
        abs_fft = np.abs(fft_power_spec) ** 2
        autocorrelation_rows = np.real(
            scifft.ifft(abs_fft, axis=0)[:freq_bins, :])  # ifft over columns

        # normalization factor
        norm_factor = np.tile(np.arange(freq_bins, 0, -1), (time_bins, 1)).T
        autocorrelation_rows = autocorrelation_rows / norm_factor

        # compute the beat spectrum
        beat_spectrum = np.mean(autocorrelation_rows, axis=1)
        # average over frequencies

        return beat_spectrum 
Example #13
Source Project: Computable   Author: ktraunmueller   File: bench_basic.py    License: MIT License 5 votes vote down vote up
def bench_random(self):
        from numpy.fft import ifft as numpy_ifft
        print()
        print('       Inverse Fast Fourier Transform')
        print('===============================================')
        print('      |     real input    |    complex input   ')
        print('-----------------------------------------------')
        print(' size |  scipy  |  numpy  |  scipy  |  numpy  ')
        print('-----------------------------------------------')
        for size,repeat in [(100,7000),(1000,2000),
                            (256,10000),
                            (512,10000),
                            (1024,1000),
                            (2048,1000),
                            (2048*2,500),
                            (2048*4,500),
                            ]:
            print('%5s' % size, end=' ')
            sys.stdout.flush()

            for x in [random([size]).astype(double),
                      random([size]).astype(cdouble)+random([size]).astype(cdouble)*1j
                      ]:
                if size > 500:
                    y = ifft(x)
                else:
                    y = direct_idft(x)
                assert_array_almost_equal(ifft(x),y)
                print('|%8.2f' % measure('ifft(x)',repeat), end=' ')
                sys.stdout.flush()

                assert_array_almost_equal(numpy_ifft(x),y)
                print('|%8.2f' % measure('numpy_ifft(x)',repeat), end=' ')
                sys.stdout.flush()

            print(' (secs for %s calls)' % (repeat))
        sys.stdout.flush() 
Example #14
Source Project: Computable   Author: ktraunmueller   File: bench_pseudo_diffs.py    License: MIT License 5 votes vote down vote up
def direct_diff(x,k=1,period=None):
    fx = fft(x)
    n = len(fx)
    if period is None:
        period = 2*pi
    w = fftfreq(n)*2j*pi/period*n
    if k < 0:
        w = 1 / w**k
        w[0] = 0.0
    else:
        w = w**k
    if n > 2000:
        w[250:n-250] = 0.0
    return ifft(w*fx).real 
Example #15
Source Project: Computable   Author: ktraunmueller   File: bench_pseudo_diffs.py    License: MIT License 5 votes vote down vote up
def direct_tilbert(x,h=1,period=None):
    fx = fft(x)
    n = len(fx)
    if period is None:
        period = 2*pi
    w = fftfreq(n)*h*2*pi/period*n
    w[0] = 1
    w = 1j/tanh(w)
    w[0] = 0j
    return ifft(w*fx) 
Example #16
Source Project: Computable   Author: ktraunmueller   File: bench_pseudo_diffs.py    License: MIT License 5 votes vote down vote up
def direct_hilbert(x):
    fx = fft(x)
    n = len(fx)
    w = fftfreq(n)*n
    w = 1j*sign(w)
    return ifft(w*fx) 
Example #17
Source Project: Computable   Author: ktraunmueller   File: bench_pseudo_diffs.py    License: MIT License 5 votes vote down vote up
def direct_shift(x,a,period=None):
    n = len(x)
    if period is None:
        k = fftfreq(n)*1j*n
    else:
        k = fftfreq(n)*2j*pi/period*n
    return ifft(fft(x)*exp(k*a)).real 
Example #18
Source Project: Computable   Author: ktraunmueller   File: test_basic.py    License: MIT License 5 votes vote down vote up
def test_definition(self):
        x = np.array([1,2,3,4+1j,1,2,3,4+2j], self.cdt)
        y = ifft(x)
        y1 = direct_idft(x)
        self.assertTrue(y.dtype == self.cdt,
                "Output dtype is %s, expected %s" % (y.dtype, self.cdt))
        assert_array_almost_equal(y,y1)

        x = np.array([1,2,3,4+0j,5], self.cdt)
        assert_array_almost_equal(ifft(x),direct_idft(x)) 
Example #19
Source Project: Computable   Author: ktraunmueller   File: test_basic.py    License: MIT License 5 votes vote down vote up
def test_definition_real(self):
        x = np.array([1,2,3,4,1,2,3,4], self.rdt)
        y = ifft(x)
        self.assertTrue(y.dtype == self.cdt,
                "Output dtype is %s, expected %s" % (y.dtype, self.cdt))
        y1 = direct_idft(x)
        assert_array_almost_equal(y,y1)

        x = np.array([1,2,3,4,5], dtype=self.rdt)
        self.assertTrue(y.dtype == self.cdt,
                "Output dtype is %s, expected %s" % (y.dtype, self.cdt))
        assert_array_almost_equal(ifft(x),direct_idft(x)) 
Example #20
Source Project: Computable   Author: ktraunmueller   File: test_basic.py    License: MIT License 5 votes vote down vote up
def test_djbfft(self):
        for i in range(2,14):
            n = 2**i
            x = list(range(n))
            y = fftpack.zfft(x,direction=-1)
            y2 = numpy.fft.ifft(x)
            assert_array_almost_equal(y,y2)
            y = fftpack.zrfft(x,direction=-1)
            assert_array_almost_equal(y,y2) 
Example #21
Source Project: Computable   Author: ktraunmueller   File: test_basic.py    License: MIT License 5 votes vote down vote up
def test_random_complex(self):
        for size in [1,51,111,100,200,64,128,256,1024]:
            x = random([size]).astype(self.cdt)
            x = random([size]).astype(self.cdt) + 1j*x
            y1 = ifft(fft(x))
            y2 = fft(ifft(x))
            self.assertTrue(y1.dtype == self.cdt,
                    "Output dtype is %s, expected %s" % (y1.dtype, self.cdt))
            self.assertTrue(y2.dtype == self.cdt,
                    "Output dtype is %s, expected %s" % (y2.dtype, self.cdt))
            assert_array_almost_equal(y1, x)
            assert_array_almost_equal(y2, x) 
Example #22
Source Project: Computable   Author: ktraunmueller   File: test_basic.py    License: MIT License 5 votes vote down vote up
def test_random_real(self):
        for size in [1,51,111,100,200,64,128,256,1024]:
            x = random([size]).astype(self.rdt)
            y1 = ifft(fft(x))
            y2 = fft(ifft(x))
            self.assertTrue(y1.dtype == self.cdt,
                    "Output dtype is %s, expected %s" % (y1.dtype, self.cdt))
            self.assertTrue(y2.dtype == self.cdt,
                    "Output dtype is %s, expected %s" % (y2.dtype, self.cdt))
            assert_array_almost_equal(y1, x)
            assert_array_almost_equal(y2, x) 
Example #23
Source Project: Computable   Author: ktraunmueller   File: test_basic.py    License: MIT License 5 votes vote down vote up
def test_djbfft(self):
        from numpy.fft import ifft as numpy_ifft
        for i in range(2,14):
            n = 2**i
            x = list(range(n))
            x1 = zeros((n,),dtype=cdouble)
            x1[0] = x[0]
            for k in range(1, n//2):
                x1[k] = x[2*k-1]+1j*x[2*k]
                x1[n-k] = x[2*k-1]-1j*x[2*k]
            x1[n//2] = x[-1]
            y1 = numpy_ifft(x1)
            y = fftpack.drfft(x,direction=-1)
            assert_array_almost_equal(y,y1) 
Example #24
Source Project: Computable   Author: ktraunmueller   File: test_basic.py    License: MIT License 5 votes vote down vote up
def test_complex(self):
        if np.dtype(np.longcomplex).itemsize == np.dtype(np.complex).itemsize:
            # longdouble == double; so fft is supported
            return

        x = np.random.randn(10).astype(np.longdouble) + \
                1j * np.random.randn(10).astype(np.longdouble)

        for f in [fft, ifft]:
            try:
                f(x)
                raise AssertionError("Type %r not supported but does not fail" %
                                     np.longcomplex)
            except ValueError:
                pass 
Example #25
Source Project: Computable   Author: ktraunmueller   File: test_basic.py    License: MIT License 5 votes vote down vote up
def test_real(self):
        if np.dtype(np.longdouble).itemsize == np.dtype(np.double).itemsize:
            # longdouble == double; so fft is supported
            return

        x = np.random.randn(10).astype(np.longcomplex)

        for f in [fft, ifft]:
            try:
                f(x)
                raise AssertionError("Type %r not supported but does not fail" %
                                     np.longcomplex)
            except ValueError:
                pass 
Example #26
Source Project: Computable   Author: ktraunmueller   File: test_basic.py    License: MIT License 5 votes vote down vote up
def test_ifft(self):
        overwritable = (np.complex128, np.complex64)
        for dtype in self.dtypes:
            self._check_1d(ifft, dtype, (16,), -1, overwritable)
            self._check_1d(ifft, dtype, (16, 2), 0, overwritable)
            self._check_1d(ifft, dtype, (2, 16), 1, overwritable) 
Example #27
Source Project: tools   Author: kastnerkyle   File: audio_tools.py    License: BSD 3-Clause "New" or "Revised" License 5 votes vote down vote up
def harvest_get_f0_candidate_from_raw_event(boundary_f0,
        fs, y_spectrum, y_length, temporal_positions, f0_floor,
        f0_ceil):
    filter_length_half = int(np.round(fs / boundary_f0 * 2))
    band_pass_filter_base = harvest_nuttall(filter_length_half * 2 + 1)
    shifter = np.cos(2 * np.pi * boundary_f0 * np.arange(-filter_length_half, filter_length_half + 1) / float(fs))
    band_pass_filter = band_pass_filter_base * shifter

    index_bias = filter_length_half
    # possible numerical issues if 32 bit
    spectrum_low_pass_filter = np.fft.fft(band_pass_filter.astype("float64"), len(y_spectrum))
    filtered_signal = np.real(np.fft.ifft(spectrum_low_pass_filter * y_spectrum))
    index_bias = filter_length_half + 1
    filtered_signal = filtered_signal[index_bias + np.arange(y_length).astype("int32")]
    negative_zero_cross = harvest_zero_crossing_engine(filtered_signal, fs)
    positive_zero_cross = harvest_zero_crossing_engine(-filtered_signal, fs)
    d_filtered_signal = filtered_signal[1:] - filtered_signal[:-1]
    peak = harvest_zero_crossing_engine(d_filtered_signal, fs)
    dip = harvest_zero_crossing_engine(-d_filtered_signal, fs)
    f0_candidate = harvest_get_f0_candidate_contour(negative_zero_cross,
            positive_zero_cross, peak, dip, temporal_positions)
    f0_candidate[f0_candidate > (boundary_f0 * 1.1)] = 0.
    f0_candidate[f0_candidate < (boundary_f0 * .9)] = 0.
    f0_candidate[f0_candidate > f0_ceil] = 0.
    f0_candidate[f0_candidate < f0_floor] = 0.
    return f0_candidate 
Example #28
Source Project: tools   Author: kastnerkyle   File: audio_tools.py    License: BSD 3-Clause "New" or "Revised" License 5 votes vote down vote up
def cheaptrick_smoothing_with_recovery(smoothed_spectrum, f0, fs, fft_size, q1):
    quefrency_axis = np.arange(fft_size) / float(fs)
    # 0 is NaN
    smoothing_lifter = np.sin(np.pi * f0 * quefrency_axis) / (np.pi * f0 * quefrency_axis)
    p = smoothing_lifter[1:int(fft_size / 2)][::-1].copy()
    smoothing_lifter[int(fft_size / 2) + 1:] = p
    smoothing_lifter[0] = 1.
    compensation_lifter = (1 - 2. * q1) + 2. * q1 * np.cos(2 * np.pi * quefrency_axis * f0)
    p = compensation_lifter[1:int(fft_size / 2)][::-1].copy()
    compensation_lifter[int(fft_size / 2) + 1:] = p
    tandem_cepstrum = np.fft.fft(np.log(smoothed_spectrum))
    tmp_spectral_envelope = np.exp(np.real(np.fft.ifft(tandem_cepstrum * smoothing_lifter * compensation_lifter)))
    spectral_envelope = tmp_spectral_envelope[:int(fft_size / 2) + 1]
    return spectral_envelope 
Example #29
Source Project: autosynch   Author: SwagLyrics   File: signal_transforms.py    License: GNU General Public License v3.0 5 votes vote down vote up
def _i_dft(magnitude_spect, phase, window_size):
    """Discrete Fourier Transformation(Synthesis) of a given spectral analysis
    via the :func:`scipy.fftpack.ifft` inverse FFT function.

    :param magnitude_spect: Magnitude spectrum.
    :type magnitude_spect: numpy.core.multiarray.ndarray
    :param phase: Phase spectrum.
    :type phase: numpy.core.multiarray.ndarray
    :param window_size: Synthesis window size.
    :type window_size: int
    :return: Time-domain signal.
    :rtype: numpy.core.multiarray.ndarray
    """
    # Get FFT Size
    fft_size = magnitude_spect.size
    fft_points = (fft_size - 1) * 2

    # Half of window size parameters
    hw_1 = int(np.floor((window_size + 1) / 2))
    hw_2 = int(np.floor(window_size / 2))

    # Initialise output spectrum with zeros
    tmp_spect = np.zeros(fft_points, dtype=complex)
    # Initialise output array with zeros
    time_domain_signal = np.zeros(window_size)

    # Compute complex spectrum(both sides) in two steps
    tmp_spect[0:fft_size] = magnitude_spect * np.exp(1j * phase)
    tmp_spect[fft_size:] = magnitude_spect[-2:0:-1] * np.exp(-1j * phase[-2:0:-1])

    # Perform the iDFT
    fft_buf = np.real(fftpack.ifft(tmp_spect))

    # Roll-back the zero-phase windowing technique
    time_domain_signal[:hw_2] = fft_buf[-hw_2:]
    time_domain_signal[hw_2:] = fft_buf[:hw_1]

    return time_domain_signal

# EOF 
Example #30
Source Project: GraphicDesignPatternByPython   Author: Relph1119   File: test_basic.py    License: MIT License 5 votes vote down vote up
def direct_idftn(x):
    x = asarray(x)
    for axis in range(len(x.shape)):
        x = ifft(x,axis=axis)
    return x