covariance.py 文件源码

def cospectrum(X, window=128, overlap=0.75, fmin=None, fmax=None, fs=None,
               phase_correction=False):
    Ne, Ns = X.shape
    number_freqs = int(window / 2)

    step = int((1.0 - overlap) * window)
    step = max(1, step)

    number_windows = (Ns - window) / step + 1
    # pre-allocation of memory
    fdata = numpy.zeros((number_windows, Ne, number_freqs), dtype=complex)
    win = numpy.hanning(window)

    # Loop on all frequencies
    for window_ix in range(int(number_windows)):

        # time markers to select the data
        # marker of the beginning of the time window
        t1 = int(window_ix * step)
        # marker of the end of the time window
        t2 = int(t1 + window)
        # select current window and apodize it
        cdata = X[:, t1:t2] * win

        # FFT calculation
        fdata[window_ix, :, :] = numpy.fft.fft(
            cdata, n=window, axis=1)[:, 0:number_freqs]

        # if(phase_correction):
        # fdata = fdata.*(exp(-sqrt(-1)*t1*( numpy.range(window)
        # ).T/window*2*pi)*numpy.ones((1,Ne))

    # Adjust Frequency range to specified range (in case it is a parameter)
    if fmin is not None:
        f = numpy.arange(0, 1, 1.0 / number_freqs) * (fs / 2.0)
        Fix = (f >= fmin) & (f <= fmax)
        fdata = fdata[:, :, Fix]

    # fdata = fdata.real
    Nf = fdata.shape[2]
    S = numpy.zeros((Ne, Ne, Nf), dtype=complex)

    for i in range(Nf):
        S[:, :, i] = numpy.dot(
            fdata[:, :, i].conj().T, fdata[:, :, i]) / number_windows

    return S