Original author: David Banas capn.freako@gmail.com
Original date: November 30, 2017
Copyright © 2017 David Banas; all rights reserved World wide.
This Jupyter notebook can be used to sanity check a batch of 4-port single-ended Touchstone files, and generate their equivalent differential step responses, for use w/ PyBERT.
%matplotlib inline
from matplotlib import pyplot as plt
import skrf as rf
import numpy as np
import scipy.signal as sig
Load the channels of interest and sanity check them.
plt.figure(figsize=(16, 8))
def sdd_21(ntwk):
"""
Given a 4-port single-ended network, return its differential throughput."
Note that the following connectivity is assumed:
- 1 ==> 2
- 3 ==> 4
"""
return 0.5*(ntwk.s21 - ntwk.s23 + ntwk.s43 - ntwk.s41)
# Load the channels of interest.
chnls = []
chnls.append(('ch1', rf.Network('../channels/802.3bj_COM_Cisco/kochuparambil_3bj_02_0913/Beth_longSmooth_THRU.s4p')))
chnls.append(('ch2', rf.Network('../channels/802.3bj_COM_Cisco/kochuparambil_3bj_02_0913/Beth_shortReflective_THRU.s4p')))
chnls.append(('ch5', rf.Network('../channels/shanbhag_01_0511/TEC_Whisper27in_THRU_G14G15.s4p')))
# Create diagonal mask, for checking passivity, below.
n = chnls[0][1]
m = np.zeros(n.s[0].shape, dtype=bool)
np.fill_diagonal(m, True)
# Check them against several criteria.
passivities = []
clrs = ['r','g','b','c','m','y']
for (lbl, ntwk), clr in zip(chnls, clrs):
passivity = map(lambda x: max(x[m]), ntwk.passivity)
# passivities.append(ntwk.is_passive()) # Apparently, not implemented, yet.
if(max(passivity) <= 1.0):
passivities.append(True)
else:
passivities.append(False)
plt.subplot(121)
plt.plot(ntwk.f / 1e9, passivity, clr, label=lbl)
plt.subplot(122)
plt.plot(ntwk.f / 1e9, sdd_21(ntwk).group_delay.flatten(), clr, label=lbl)
#Plot passivities.
plt.subplot(121)
plt.title("Passivity Plots")
plt.xlabel("f (GHz)")
plt.grid()
plt.legend(loc='upper right')
#Plot group delays.
plt.subplot(122)
plt.title("Group Delay")
plt.xlabel("f (GHz)")
plt.grid()
plt.legend(loc='upper right')
# Print a summary of the results.
print("{:10s} {:10s}".format('Channel','Passive'))
print('_'*21)
for ((lbl, ntwk), passive) in zip(chnls, passivities):
print("{:^10s} {:^10s}".format(lbl, str(passive)))
All channels appear passive. And none have any serious anomalies in their group delay plots.
Convert to mixed mode and display Sdd[2,1].
# Plot their differential insertion losses.
sdd_21s = []
plt.figure(figsize=(16, 8))
for (lbl, ntwk), clr in zip(chnls, clrs):
H = sdd_21(ntwk)
sdd_21s.append(H)
plt.semilogx(ntwk.f, 20 * np.log10(abs(H.s[:,0,0])), clr, label=lbl)
plt.title("SDD[2,1]")
plt.xlabel("f (Hz)")
plt.ylabel("|Sdd21| (dB)")
plt.grid()
plt.legend(loc='upper right')
plt.axis(ymin=-40)
No real surprises, here. So far, so good.
# Here, I calculate the impulse response myself,
# as opposed to using the 'Network.s_time()' function provided by skrf,
# because I've found I get better results.
ts = []
fs = []
hs = []
ss = []
plt.figure(figsize=(16, 8))
for (lbl, _), sdd_21, clr in zip(chnls, sdd_21s, clrs):
# Form frequency vector.
f = sdd_21.f
fmin = f[0]
if(fmin == 0): # Correct, if d.c. point was included in original data.
fmin = f[1]
fmax = f[-1]
f = np.arange(fmin, fmax + fmin, fmin)
F = rf.Frequency.from_f(f / 1e9) ## skrf.Frequency.from_f() expects its argument to be in units of GHz.
print("{} frequencies: {}".format(lbl, F))
# Form impulse response from frequency response.
sdd_21 = sdd_21.interpolate_from_f(F)
H = sdd_21.s[:,0,0]
H = np.concatenate((H, np.conj(np.flip(H[:-1], 0)))) # Forming the vector that fft() would've outputted.
H = np.pad(H, (1,0), 'constant', constant_values=1.0) # Presume d.c. value = 1.
h = np.real(np.fft.ifft(H))
h /= np.abs(h.sum()) # Equivalent to assuming that step response settles at 1.
# Form step response from impulse response.
s = np.cumsum(h)
# Form time vector.
t0 = 1. / (2. * fmax) # Sampling interval = 1 / (2 fNyquist).
t = np.array([n * t0 for n in range(len(h))])
# Save results.
ts.append(t)
fs.append(f)
hs.append(h)
ss.append(s)
# Plot results.
plt.plot(t * 1e9, h / t0 * 1e-9, clr, label=lbl)
# Annotate the plot.
plt.title("Impulse Response")
plt.xlabel("t (ns)")
plt.ylabel("h(t) (V/ns)")
plt.grid()
plt.legend(loc='upper right')
plt.axis(xmin=0, xmax=10);
Note that the locations of the peaks in the plot, above, agree with the group delays plotted, further above. This lends confidence to the code used to extract these impulse responses from the S-parameter data.
Further, note that the magnitudes are all very reasonable.
Channel, ch2, has some troubling artifacts just before the initial attack. This is, probably, due to insufficient bandwidth and an insufficient number of frequency sampling points for a channel this short and reflective.
plt.figure(figsize=(16, 8))
for (lbl, _), t, s, clr in zip(chnls, ts, ss, clrs):
plt.plot(t * 1e9, s, clr, label=lbl)
plt.title("Step Response")
plt.xlabel("t (ns)")
plt.ylabel("s(t) (V)")
plt.grid()
plt.legend(loc='upper left');
# plt.axis(xmin=1.8, xmax=2.0)
Channels, ch1 and ch2, are clearly problematic. ch5 is the cleanest.
# Save the step responses.
for (lbl, _), t, s in zip(chnls, ts, ss):
with open(lbl+'_s.csv', 'wt') as file:
for x, y in zip(t, s):
print >> file, "{:014.12f}, {:05.3f}".format(x, y)
# Take them back into the frequency domain and compare to original data.
hps = []
plt.figure(figsize=(16, 8))
for f, s, (lbl, _), Href, clr in zip(fs, ss, chnls, sdd_21s, clrs):
h = np.diff(s)
h = np.pad(h, (1, 0), 'constant', constant_values=0)
H = np.fft.fft(h)
# Using [1], as opposed to [0], to accomodate ch1. (See, below.)
# I'm assuming the strange behavior in ch1 is due to me forcing d.c. values of 1.
H *= abs(Href.s[:,0,0][1]) / abs(H[1]) # Normalize the "d.c." levels.
plt.semilogx(Href.f, 20 * np.log10(abs(Href.s[:,0,0])), clr+'--', label=lbl+'_ref')
plt.semilogx(f, 20 * np.log10(abs(H[:len(f)])), clr, label=lbl)
plt.title("SDD[2,1]")
plt.xlabel("f (Hz)")
plt.ylabel("|Sdd21| (dB)")
plt.grid()
plt.legend(loc='upper right')
plt.axis(ymin=-40);
With the single exception of ch1, the comparisons look great. I claim we have faithful transformations from S-parameter data to step responses.
Wrt/ ch1, note that the only anomaly occurs in the very first frequency point (50 MHz) and probably won't cause any issues in PyBERT. (It's the high frequency content that matters most, when trying to open eyes in PyBERT.)
Note: Don't forget that both ch1 and ch2 have exhibited signs of a significant lack of numerical fidelity throughout this analysis.
(The ch1 / ch2 results are a bit of a surprise, actually, because I pulled them from the 802.3bj representative channel library. I wonder if those folks are aware of this.)