Merge 93df9d4 into e2f9497

examples/AWG8_examples/bounce_correction_test.py

@@ -0,0 +1,107 @@
+
+import numpy as np
+import matplotlib as mpl
+import matplotlib.pyplot as plt
+from scipy import signal
+
+
+import pycqed.measurement.kernel_functions_ZI as ZI_kern
+
+
+mpl.rcParams['font.size'] = 12
+mpl.rcParams['legend.fontsize'] = 12
+mpl.rcParams['figure.titlesize'] = 'medium'
+
+
+# Settings
+
+
+fs = 2.4e9
+time_start = -50e-9
+time_start = np.around(time_start*fs)/fs
+time_end = 50e-9
+time = np.arange(time_start, time_end, 1/fs)
+
+delay = 10.1e-9
+amplitude = 0.1
+
+
+# Construct impulse_response
+impulse = np.zeros(len(time))
+zero_ind = np.argmin(np.abs(time))
+impulse[zero_ind] = 1.0
+delay_ind = np.argmin(np.abs(time-delay))
+impulse_response = np.copy(impulse)
+impulse_response[delay_ind] = amplitude
+
+
+# Derive step response
+step = np.zeros(len(time))
+step[time >= 0.0] = 1.0
+step_response = signal.lfilter(impulse_response[zero_ind:], 1.0, step)
+
+
+# Compute ideal inverted filter kernel
+a = ZI_kern.ideal_inverted_fir_kernel(impulse_response, zero_ind)
+
+
+# Apply ideal inverted filter to impulse response and step response
+impulse_response_corr = signal.lfilter(a, 1.0, impulse_response)
+step_response_corr = signal.lfilter(a, 1.0, step_response)
+
+# Apply hardware-friendly filter to impulse response and step response
+impulse_response_corr_hw = ZI_kern.multipath_bounce_correction(impulse_response, round(delay*fs), -amplitude)
+step_response_corr_hw = ZI_kern.multipath_bounce_correction(step_response, round(delay*fs), -amplitude)
+
+
+# Plot impulse response comparison
+plt.figure(1, figsize=(7,10))
+
+plt.subplot(3, 1, 1)
+plt.plot(time*1e9, impulse_response)
+plt.xlabel('Time, t (ns)')
+plt.ylabel('Amplitude (a.u)')
+plt.title('(a) Impulse response')
+
+plt.subplot(3, 1, 2)
+plt.plot(time*1e9, impulse_response_corr)
+plt.xlabel('Time, t (ns)')
+plt.ylabel('Amplitude (a.u)')
+plt.title('(b) Ideal corrected impulse response')
+
+plt.subplot(3, 1, 3)
+plt.plot(time*1e9, impulse_response_corr_hw)
+plt.xlabel('Time, t (ns)')
+plt.ylabel('Amplitude (a.u)')
+plt.title('(c) Harware-corrected impulse response')
+
+plt.tight_layout()
+plt.savefig('impulse_response.png',dpi=600,bbox_inches='tight')
+plt.show()
+
+
+# Plot step response comparison
+plt.figure(1, figsize=(7,10))
+
+plt.subplot(3, 1, 1)
+plt.plot(time*1e9, step_response)
+plt.xlabel('Time, t (ns)')
+plt.ylabel('Amplitude (a.u)')
+plt.title('(a) Step response')
+
+plt.subplot(3, 1, 2)
+plt.plot(time*1e9, step_response_corr)
+plt.xlabel('Time, t (ns)')
+plt.ylabel('Amplitude (a.u)')
+plt.title('(b) Ideal corrected step response')
+
+plt.subplot(3, 1, 3)
+plt.plot(time*1e9, step_response_corr_hw)
+plt.xlabel('Time, t (ns)')
+plt.ylabel('Amplitude (a.u)')
+plt.title('(c) Harware-corrected step response')
+
+plt.tight_layout()
+plt.savefig('step_response.png',dpi=600,bbox_inches='tight')
+plt.show()
+

pycqed/measurement/kernel_functions_ZI.py

@@ -168,3 +168,37 @@ def multipath_filter2(sig, alpha, k, paths):
             duf = np.append(duf, tpl * acc[j])
     duf = duf[0:sig.size]
     return sig + k * (duf - sig)
+
+def multipath_bounce_correction(sig, delay, amp, paths = 8, bufsize = 128):
+    """
+    This function simulates a possible FPGA implementation of a first-order bounce correction filter.
+    The signal (sig) is assumed to be a numpy array representing a wavefomr with sampling rate 2.4 GSa/s.
+    The delay is specified in number of samples. It needs to be an interger.
+    The amplitude (amp) of the bounce is specified relative to the amplitude of the input signal.
+    It is constrained to be smaller than 1. The amplitude is represented as a 18-bit fixed point number on the FPGA.
+    """
+    assert 0 <= delay < bufsize-8, "The maximulm delay is limitted to 120 (bufsize-8) samples to save hardware resources."
+    assert -1 <= amp < 1, "The amplitude needs to be between -1 and 1."
+
+    sigout = np.zeros(len(sig))
+    buffer = np.zeros(bufsize)
+
+    amp_hw = coef_round(amp)
+
+    # iterate in steps of eight samples through the input signal to simulate the implementation with parallel paths on the FPGA
+    for i in range(0, len(sig), paths):
+        buffer[paths:] = buffer[:-paths]
+        buffer[:paths] = sig[i:i+8]
+        sigout[i:i+8] = sig[i:i+8] + amp_hw*buffer[delay:delay+8]
+    return sigout
+
+
+def ideal_inverted_fir_kernel(impulse_response, zero_ind=0):
+    """
+    This function computes the ideal inverted FIR filter kernel for a given impulse_response.
+    The argument zero_ind provides the index corresponding to time t=0 within the impulse_response array.
+    """
+    f = np.fft.fft(impulse_response)
+    impulse_response_inv = np.fft.ifft(1/f)
+    impulse_response_inv_re_trunc= np.real(impulse_response_inv)[zero_ind:]
+    return impulse_response_inv_re_trunc