Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
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
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)
# 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)
# 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]
_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)
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()
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)
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)