Skip to content
Snippets Groups Projects
Commit c249d74b authored by Kipp Cannon's avatar Kipp Cannon
Browse files

linear_phase_fir_kernel_to_minimum_phase_whitening_fir_kernel():

- remove superfluous array slice operations
- future-proof integer divisions
- "inline" some algebra to reduce swig interaction overhead
parent 07b5f570
No related branches found
No related tags found
No related merge requests found
......@@ -533,15 +533,14 @@ class PSDFirKernel(object):
)
lal.COMPLEX16TimeFreqFFT(absX, kernel_tseries, self.fwdplan)
absX.data.data[:] = numpy.roll(abs(absX.data.data), -working_length / 2 + 1)[:] * sample_rate
absX.data.data[:] = numpy.roll(abs(absX.data.data), -working_length // 2 + 1) * sample_rate
#
# compute the cepstrum of the kernel
# (i.e., the iFFT of the log of the abs of the FFT of the kernel)
#
logabsX.data.data = numpy.log(absX.data.data)
logabsX.data.data[:] = numpy.roll(logabsX.data.data, +working_length / 2 + 1)[:]
logabsX.data.data[:] = numpy.roll(numpy.log(absX.data.data), +working_length // 2 + 1)
lal.COMPLEX16FreqTimeFFT(cepstrum, logabsX, self.revplan)
cepstrum.data.data /= sample_rate
......@@ -551,26 +550,25 @@ class PSDFirKernel(object):
sgn = scipy.ones(working_length)
sgn[0] = 0.
sgn[(len(sgn)+1)/2] = 0.
sgn[(len(sgn)+1)/2:] *= -1
sgn[(len(sgn) + 1) // 2] = 0.
sgn[(len(sgn) + 1) // 2:] *= -1
cepstrum.data.data *= sgn
#
# compute theta
#
lal.COMPLEX16TimeFreqFFT(theta, cepstrum, self.fwdplan)
theta.data.data *= -1.j
theta.data.data[:] = numpy.roll(theta.data.data, -working_length / 2)[:] * sample_rate
theta.data.data[:] = numpy.roll(-1.j * theta.data.data, -working_length // 2) * sample_rate
#
# compute minimum phase kernel
#
absX.data.data[:] = absX.data.data * scipy.exp(1.j * theta.data.data)[:]
absX.data.data[:] = numpy.roll(absX.data.data, +working_length / 2 + 1)[:] / sample_rate
absX.data.data[:] = numpy.roll(absX.data.data * scipy.exp(1.j * theta.data.data), +working_length // 2 + 1) / sample_rate
lal.COMPLEX16FreqTimeFFT(min_phase_kernel, absX, self.revplan)
kernel = numpy.real(min_phase_kernel.data.data)
kernel = min_phase_kernel.data.data.real
#
# this kernel needs to be reversed to follow conventions used with the
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment