How to use the pycbc.fft.ifft function in PyCBC

To help you get started, we’ve selected a few PyCBC examples, based on popular ways it is used in public projects.

Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.

github gwastro / pycbc / pycbc / psd / estimate.py View on Github external
Returns
    -------
    interpolated series : FrequencySeries
        A new FrequencySeries that has been interpolated.
    """
    series = FrequencySeries(series, dtype=complex_same_precision_as(series), delta_f=series.delta_f)

    N = (len(series) - 1) * 2
    delta_t = 1.0 / series.delta_f / N

    new_N = int(1.0 / (delta_t * delta_f))
    new_n = new_N // 2 + 1

    series_in_time = TimeSeries(zeros(N), dtype=real_same_precision_as(series), delta_t=delta_t)
    ifft(series, series_in_time)

    padded_series_in_time = TimeSeries(zeros(new_N), dtype=series_in_time.dtype, delta_t=delta_t)
    padded_series_in_time[0:N//2] = series_in_time[0:N//2]
    padded_series_in_time[new_N-N//2:new_N] = series_in_time[N//2:N]

    interpolated_series = FrequencySeries(zeros(new_n), dtype=series.dtype, delta_f=delta_f)
    fft(padded_series_in_time, interpolated_series)

    return interpolated_series
github gwastro / pycbc / pycbc / filter / qtransform.py View on Github external
start = int((f0 - (f0 / qprime)) * fseries.duration)
    end = int(start + window_size)
    center = (start + end) // 2

    windowed = fseries[start:end] * (1 - xfrequencies ** 2) ** 2 * norm

    tlen = (len(fseries)-1) * 2
    windowed.resize(tlen)
    windowed.roll(-center)

    # calculate the time series for this q -value
    windowed = FrequencySeries(windowed, delta_f=fseries.delta_f,
                            epoch=fseries.start_time)
    ctseries = TimeSeries(zeros(tlen, dtype=numpy.complex128),
                            delta_t=fseries.delta_t)
    ifft(windowed, ctseries)

    if return_complex:
        return ctseries
    else:
        energy = ctseries.squared_norm()
        medianenergy = numpy.median(energy.numpy())
        return  energy / float(medianenergy)
github gwastro / pycbc / test / fft_base.py View on Github external
# Now the same for FFT(IFFT(random))
    if dtype(outarr).kind == 'c':
        outarr._data[:] = randn(len(outarr))+1j*randn(len(outarr))
        # If we're going to do a HC2R transform we must worry about DC/Nyquist imaginary
        if dtype(inarr).kind == 'f':
            outarr._data[0] = real(outarr[0])
            if (len(inarr)%2)==0:
                outarr._data[len(outarr)-1] = real(outarr[len(outarr)-1])
    else:
        outarr._data[:] = randn(len(outarr))
    inarr.clear()
    outcopy = type(outarr)(outarr)
    if type(outarr) == pycbc.types.Array:
        outcopy *= len(inarr)
    with tc.context:
        pycbc.fft.ifft(outarr, inarr)
        pycbc.fft.fft(inarr, outarr)
        emsg="FFT(IFFT(random)) did not reproduce original array to within tolerance {0}".format(tol)
        if isinstance(outcopy,ts) or isinstance(outcopy,fs):
            tc.assertTrue(outcopy.almost_equal_norm(outarr,tol=tol,dtol=tol),
                          msg=emsg)
        else:
            tc.assertTrue(outcopy.almost_equal_norm(outarr,tol=tol),
                          msg=emsg)
github gwastro / pycbc / test / fft_base.py View on Github external
# the whole point is to vary the input/output in some way that should cause
    # an exception.
    if other_args is None:
        other_args = {}
    tc = test_case
    with tc.context:
        outty = type(outarr)
        outzer = pycbc.types.zeros(len(outarr))
        # If we give an output array that is wrong only in length, raise ValueError:
        out_badlen = outty(pycbc.types.zeros(len(outarr)+1),dtype=outarr.dtype,**other_args)
        args = [inarr, out_badlen]
        tc.assertRaises(ValueError,pycbc.fft.ifft,*args)
        # If we give an output array that has the wrong precision, raise ValueError:
        out_badprec = outty(outzer,dtype=_other_prec[dtype(outarr).type],**other_args)
        args = [inarr,out_badprec]
        tc.assertRaises(ValueError,pycbc.fft.ifft,*args)
        # If we give an output array that has the wrong kind (real or complex) but
        # correct precision, then raise a ValueError.  Here we must adjust the kind
        # of the *input* array, not output.  But that makes it hard, because the 'other_args'
        # parameter will be wrong for that.  Very hacky, but oh well...
        new_args = other_args.copy()
        if new_args != {}:
            try:
                delta = new_args.pop('delta_t')
                new_args.update({'delta_f' : delta})
            except KeyError:
                delta = new_args.pop('delta_f')
                new_args.update({'delta_t' : delta})
        in_badkind = type(inarr)(pycbc.types.zeros(len(inarr)),
                                 dtype=_bad_dtype[dtype(outarr).type],
                                 **new_args)
        args = [in_badkind, outarr]
github gwastro / pycbc / pycbc / filter / matchedfilter.py View on Github external
_q = out
    else:
        raise TypeError('Invalid Output Vector: wrong length or dtype')

    correlate(htilde[kmin:kmax], stilde[kmin:kmax], qtilde[kmin:kmax])

    if psd is not None:
        if isinstance(psd, FrequencySeries):
            if psd.delta_f == stilde.delta_f :
                qtilde[kmin:kmax] /= psd[kmin:kmax]
            else:
                raise TypeError("PSD delta_f does not match data")
        else:
            raise TypeError("PSD must be a FrequencySeries")

    ifft(qtilde, _q)

    if h_norm is None:
        h_norm = sigmasq(htilde, psd, low_frequency_cutoff, high_frequency_cutoff)

    norm = (4.0 * stilde.delta_f) / sqrt( h_norm)

    return (TimeSeries(_q, epoch=stilde._epoch, delta_t=stilde.delta_t, copy=False),
           FrequencySeries(qtilde, epoch=stilde._epoch, delta_f=stilde.delta_f, copy=False),
           norm)
github gwastro / pycbc / examples / waveform / plot_fd_td.py View on Github external
import pylab
from pycbc import types, fft, waveform

# Get a time domain waveform
hp, hc = waveform.get_td_waveform(approximant="EOBNRv2",
                             mass1=6, mass2=6, delta_t=1.0/4096, f_lower=40)

# Get a frequency domain waveform
sptilde, sctilde = waveform. get_fd_waveform(approximant="TaylorF2",
                             mass1=6, mass2=6, delta_f=1.0/4, f_lower=40)

# FFT it to the time-domain
tlen = int(1.0 / hp.delta_t / sptilde.delta_f)
sptilde.resize(tlen/2 + 1)
sp = types.TimeSeries(types.zeros(tlen), delta_t=hp.delta_t)
fft.ifft(sptilde, sp)

pylab.plot(sp.sample_times, sp, label="TaylorF2 (IFFT)")
pylab.plot(hp.sample_times, hp, label="EOBNRv2")

pylab.ylabel('Strain')
pylab.xlabel('Time (s)')
pylab.legend()
pylab.show()
github gwastro / pycbc / pycbc / filter / resample.py View on Github external
timeseries = Array(timeseries, copy=False)

        flen = len(cseries) / 2 + 1
        ftype = complex_same_precision_as(timeseries)

        cfreq = zeros(flen, dtype=ftype)
        tfreq = zeros(flen, dtype=ftype)

        fft(Array(cseries), cfreq)
        fft(Array(timeseries), tfreq)

        cout = zeros(flen, ftype)
        out = zeros(len(timeseries), dtype=timeseries)

        correlate(cfreq, tfreq, cout)
        ifft(cout, out)

        return out.numpy()  / len(out)
github gwastro / pycbc / pycbc / vetoes / powert_chisq.py View on Github external
def __init__ (self, num_bins, segs, psd, fmin):
        self.nbins = num_bins
        self.kmin = int(fmin / psd.delta_f)
        self.asd = psd ** 0.5  
        self.segs = []
        
        self.N = (len(psd)-1)*2
        self.dt = 1.0 / (self.N * psd.delta_f)
         
        # Pregenerate the whitened h(t) segments 
        for stilde in segs:
            ht = TimeSeries(zeros(self.N), delta_t=self.dt, dtype=self.asd.dtype)    
            stilde = stilde / self.asd
            stilde[0:self.kmin].clear()  
            ifft(stilde, ht)
            self.segs.append(ht)