# Time domain reflectometry, measurement vs simulation¶

This example demonstrate the use of frequency to time domain transformation by comparing measurements and simulations of a microstripline and a microstripline with stepped impedance sections.

The simulation data is generated by skrf using a simple transmission line model for connectors and each impedance section. To achieve a reasonnable agreement between measured and simulated data, the dielectric permittivity as well as the connector impedance and delay are extracted by optimization. The code for the simulation and parameters optimization is given at the end of the example.

## Data preparation¶

### Setup¶

In [1]:

%matplotlib inline
import skrf
import numpy as np
from numpy import real, log, log10, sum, absolute, pi, sqrt
from scipy.optimize import minimize
import matplotlib.pyplot as plt
from IPython.display import *
skrf.stylely()


The measurement was performed the 19th February 2018 on an Anritsu MS46524B 20GHz Vector Network Analyser. The setup is a linear frequency sweep from 1MHz to 10GHz with 1MHz step, 10000 points. If bandwidth 1kHz, 0dBm output power, no smoothing, no averaging. Full two-port calibration with eCal kit.

Considerations about time resolution and range limitation to avoid alias response are at the end of this document.

In [2]:

#load all measurement and simulation data into a dictionary


### DC point extrapolation¶

The measured and simulated data are available from 1 Mhz to 10 GHz with 1 MHz step (harmonic sweep), thus we need to extrapolate the dc point.

In [3]:

names = ['P1-MSL_Stepped_140-P2',
'P1-MSL_Thru_100-P2']

meas_dc_ext = meas.copy()
simu_dc_ext = simu.copy()
for k, n in enumerate(names):
meas_dc_ext[n] = meas_dc_ext[n].extrapolate_to_dc(kind='linear')
simu_dc_ext[n] = simu_dc_ext[n].extrapolate_to_dc(kind='linear')



## Microstripline¶

MSL_Thru_100 is a $$L1$$ long, $$W1$$ wide, $$T$$ thick copper microstripline on a $$H$$ height substrate with bottom ground plane.

Parameter Value
$$L1$$ 100mm
$$W1$$ 3.0mm
$$T$$ 50umm
$$H$$ 1.5mm
Connector Cinch, 142-0701-851
Substrate FR-4
In [4]:

Image('tdr_measurement_vs_simulation/figures/MSL_100.jpg', width='50%')

Out[4]:


### Measurement vs simulation comparison¶

In [5]:

plt.figure()
plt.subplot(2,2,1)
plt.title('Time')
meas_dc_ext['P1-MSL_Thru_100-P2'].s11.plot_z_time_step(window='hamming', label='meas')
simu_dc_ext['P1-MSL_Thru_100-P2'].s11.plot_z_time_step(window='hamming', label='simu')
plt.xlim((-2, 3))

plt.subplot(2,2,2)
plt.title('Frequency')
meas_dc_ext['P1-MSL_Thru_100-P2'].s11.plot_s_db(label='meas')
simu_dc_ext['P1-MSL_Thru_100-P2'].s11.plot_s_db(label='simu')

plt.subplot(2,2,3)
z0 = 50
ymea[ymea ==  1.] =  1. + 1e-12  # solve numerical singularity
ymea[ymea == -1.] = -1. + 1e-12  # solve numerical singularity
ymea = z0 * (1+ymea) / (1-ymea)
ysim[ysim ==  1.] =  1. + 1e-12  # solve numerical singularity
ysim[ysim == -1.] = -1. + 1e-12  # solve numerical singularity
ysim = z0 * (1+ysim) / (1-ysim)
plt.xlabel('Time (ns)')
plt.ylabel('Relative error (%)')
plt.plot(t*1e9, 100*(ysim-ymea)/ymea)
plt.xlim((-2, 3))

plt.subplot(2,2,4)
delta = simu_dc_ext['P1-MSL_Thru_100-P2'].s_db[:,0,0] - meas_dc_ext['P1-MSL_Thru_100-P2'].s_db[:,0,0]
f = simu_dc_ext['P1-MSL_Thru_100-P2'].f * 1e-9
plt.xlabel('Frequency (GHz)')
plt.ylabel('Delta (dB)')
plt.plot(f, delta)
plt.ylim((-20,20))

plt.tight_layout()
plt.show()


Surprisingly, the time domain results shows a very good agreement, within ±1%, while the frequency domain results exhibit a reasonnable agreement only on the lower half of the frequencies.

This is because the time domain base shape is mostly impacted by the low frequency data.

There is a small offset between the time domain data, sign that the DC point is a bit different. The inductive peaks caused by the connector-to-microstripline transition are clearly visible and their simulation model agree well with measurements.

To go further, it would be possible to measure up to a bigger frequency span, thus increasing the resolution and build the dut on a more frequency stable substrate, but this would be more costly. A cross section of the transmission line could also improve the knowledge of the actual geometry, including manufacturing tolerances.

## Stepped impedance microstripline¶

MSL_Stepped_100 is a stepped microstripline made of $$T$$ thick copper on a $$H$$ height substrate with bottom ground plane. The section 1, 2, 3 and 4 from the left to the right, are described with length $$Lx$$ and width $$Wx$$, where $$x$$ is the section number.

Parameter Value
$$L1$$ 50mm
$$W1$$ 3.0mm
$$L2$$ 20mm
$$W2$$ 8.0mm (capacitive)
$$L3$$ 20mm
$$W3$$ 1.0mm (inductive)
$$L4$$ 50mm
$$W4$$ 3.0mm
$$T$$ 50umm
$$H$$ 1.5mm
Connector Cinch, 142-0701-851
Substrate FR-4
In [6]:

Image('tdr_measurement_vs_simulation/figures/MSL_Stepped_140.jpg', width='75%')

Out[6]:


### Measurement vs simulation comparison¶

In [7]:

plt.figure()
plt.subplot(2,2,1)
plt.title('Time')
meas_dc_ext['P1-MSL_Stepped_140-P2'].s11.plot_z_time_step(window='hamming', label='measurement')
simu_dc_ext['P1-MSL_Stepped_140-P2'].s11.plot_z_time_step(window='hamming', label='simulation')
plt.xlim((-1, 3))

plt.subplot(2,2,2)
plt.title('Frequency')
meas_dc_ext['P1-MSL_Stepped_140-P2'].s11.plot_s_db(label='measurement')
simu_dc_ext['P1-MSL_Stepped_140-P2'].s11.plot_s_db(label='simulation')

plt.subplot(2,2,3)
z0 = 50
ymea[ymea ==  1.] =  1. + 1e-12  # solve numerical singularity
ymea[ymea == -1.] = -1. + 1e-12  # solve numerical singularity
ymea = z0 * (1+ymea) / (1-ymea)
ysim[ysim ==  1.] =  1. + 1e-12  # solve numerical singularity
ysim[ysim == -1.] = -1. + 1e-12  # solve numerical singularity
ysim = z0 * (1+ysim) / (1-ysim)
plt.xlabel('Time (ns)')
plt.ylabel('Relative error (%)')
plt.plot(t*1e9, 100*(ysim-ymea)/ymea)
plt.xlim((-2, 3))

plt.subplot(2,2,4)
delta = simu_dc_ext['P1-MSL_Stepped_140-P2'].s_db[:,0,0] - meas_dc_ext['P1-MSL_Stepped_140-P2'].s_db[:,0,0]
f = simu_dc_ext['P1-MSL_Stepped_140-P2'].f * 1e-9
plt.xlabel('Frequency (GHz)')
plt.ylabel('Delta (dB)')
plt.plot(f, delta)
plt.ylim((-10,10))

plt.tight_layout()
plt.show()


Both time domain and the frequency domain results shows a reasonnable agreement, within ±5% for time domain. The frequency domain results exhibit a good agreement, within ±dB only on the lower half of the frequencies, and the upper frequencies shows better agreement than in the previous microstripline case.

An explanation to the better frequency domain agreement is that the impedance steps induce such a discontinuities that more energi is reflected back, leading to a shorter effective length and an increased accuracy (more signal).

The capacitive and inductive section of the stepped line are clearly visible and their simulation model agree well with measurements. The first connector effect is hidden by the plot scale, while the second is masked by the reflexions.

To go further, it would again be possible to measure up to a bigger frequency span, thus increasing the resolution and build the dut on a more frequency stable substrate, but this would be more costly. A cross section of the transmission line could also improve the knowledge of the actual geometry, including manufacturing tolerances.

Eventually, the stepped discontinuities could be made smaller, to produce a smaller effect on the overall measurement (each discontinuity has an influence on the following time domain signal shape).

## Parameters optimization¶

### Dielectric effective relative permittivity and loss tangent characterisation based on multiline method¶

Only two lines with different lengths are required for dielectric permittivity and loss tangent characterisation. Since we have measured the reflects too, we will use them instead of fake reflects. We don’t have the switch terms either, but as we only extract the dielectric permittivity and loss tangent rather that doing a real calibration, there is no problem with that.

Multiline calibration algorithm get rid of connectors effects by using multiple lengths of lines. The docstring explain that At every frequency point there should be at least one line pair that has a phase difference that is not 0 degree or a multiple of 180 degree otherwise calibration equations are singular and accuracy is very poor. These conditions will not be met with the chosen lines combination, but we will still be able to get a decent estimation of the dielectric permittivity and loss tangent.

In [8]:

Image('tdr_measurement_vs_simulation/figures/MSL_100.jpg', width='50%')

Out[8]:

In [9]:

Image('tdr_measurement_vs_simulation/figures/MSL_200.jpg', width='100%')

Out[9]:

In [10]:

Image('tdr_measurement_vs_simulation/figures/MSL_Short_50.jpg', width='25%')

Out[10]:


#### Perform NISTMultilineTRL algorithm¶

Discard the No switch term provided warning.

In [11]:

# Perform NISTMultilineTRL algorithm
line100mm = meas['P1-MSL_Thru_100-P2']
line200mm = meas['P1-MSL_Thru_200-P2']
short50mm = skrf.network.two_port_reflect(meas['P1-MSL_Short_50'], meas['P2-MSL_Short_50'])
measured  = [line100mm, short50mm, line200mm]
Grefls    = [-1]
lengths   = [100e-3, 200e-3] # in meter
offset    = [50e-3] # in meter
cal = skrf.NISTMultilineTRL(measured, Grefls, lengths, er_est=4.5, refl_offset=offset)

/home/docs/checkouts/readthedocs.org/user_builds/scikit-rf/conda/stable/lib/python3.5/site-packages/scikit_rf-0.14.9-py3.5.egg/skrf/calibration/calibration.py:2016: UserWarning: No switch terms provided
warn('No switch terms provided')


#### Relative dielectric permittivity and loss tangent¶

The NISTMultilineTRL calibration give the frequency-dependant effective relative dielectric of the geometry, which is a mix between air and substrate relative permittivity.

Unfortunately, a single value of relative dielectric permittivity of the substrate at a given frequency is required for the microstripline media simulation, instead of the frequency-dependant effective value of the geometry we get from the calibration.

To overcome this difficult situation, the microstripline model which will be used is fitted by optimization on the calibration results in such a way the difference between the frequency-dependant effective relative dielectric permittivity of both dataset is minimised.

Additionnaly, a weighted contribution of the dielectric loss tangent is inserted in the optimization to minimise the difference between the measured modelled attenuation.

The optimizatiion results are dielectric $$\epsilon_r$$ and $$\tan{\delta}$$ at 1 GHz and will be used in the dielectric dispersion model.

In [12]:

# frequency axis
freq  = line100mm.frequency
f     = line100mm.frequency.f
f_ghz = line100mm.frequency.f/1e9

# the physical dimensions of the lines are known by design (neglecting manufacturing tolerances)
W   = 3.00e-3
H   = 1.51e-3
T   = 50e-6
L     = 0.1

# calibration results to compare against
ep_r_mea = cal.er_eff.real
A_mea    = 20/log(10)*cal.gamma.real

# starting values for the optimizer
A     = 0.0
f_A   = 1e9
ep_r0 = 4.5
tanD0 = 0.02
f_epr_tand  = 1e9
x0 = [ep_r0, tanD0]

# function to be minimised
def model(x, freq, ep_r_mea, A_mea, f_ep):
ep_r, tanD = x[0], x[1]
m = skrf.media.MLine(frequency=freq, Z0=50, w=W, h=H, t=T,
ep_r=ep_r, mu_r=1, rho=1.712e-8, tand=tanD, rough=0.15e-6,
f_low=1e3, f_high=1e12, f_epr_tand=f_ep,
diel='djordjevicsvensson', disp='kirschningjansen')
ep_r_mod = m.ep_reff_f.real
A_mod = m.alpha * 20/log(10)
return sum((ep_r_mod - ep_r_mea)**2)  + 0.1*sum((A_mod - A_mea)**2)

# run optimizer
res = minimize(model, x0, args=(freq, ep_r_mea, A_mea, f_epr_tand),
bounds=[(4.0, 5.0), (0.001, 0.1)])

# get the results and print the results
ep_r, tanD = res.x[0], res.x[1]
print('epr={:.3f}, tand={:.4f} at {:.1f} GHz.'.format(ep_r, tanD, f_epr_tand * 1e-9))

# build the corresponding media
m = skrf.media.MLine(frequency=freq, Z0=50, w=W, h=H, t=T,
ep_r=ep_r, mu_r=1, rho=1.712e-8, tand=tanD, rough=0.15e-6,
f_low=1e3, f_high=1e12, f_epr_tand=f_epr_tand,
diel='djordjevicsvensson', disp='kirschningjansen')

epr=4.403, tand=0.0176 at 1.0 GHz.


Calibration based values are plotted against modelled value as a sanity check.

In [13]:

plt.figure()
plt.subplot(2,2,1)
plt.xlabel('Frequency [GHz]')
plt.ylabel('$\epsilon_{r,eff}$')
plt.plot(f_ghz, ep_r_mea, label='measured')
plt.plot(f_ghz, m.ep_reff_f.real, label='model')
plt.legend()

plt.subplot(2,2,2)
plt.xlabel('Frequency [GHz]')
plt.ylabel('A (dB/m)')
plt.plot(f_ghz, A_mea, label='measured')
A_mod = 20/log(10)*m.alpha
plt.plot(f_ghz, A_mod, label='model')
plt.legend()

plt.subplot(2,2,3)
plt.xlabel('Frequency [GHz]')
plt.ylabel('$\epsilon_{r,eff}$ error [%]')
rel_err = 100 * ((ep_r_mea - m.ep_reff_f.real)/ep_r_mea)
plt.plot(f_ghz, rel_err)
plt.ylim((-2,2))

plt.subplot(2,2,4)
plt.xlabel('Frequency [GHz]')
plt.ylabel('$A$ error [%]')
rel_err = 100 * ((A_mea - A_mod)/A_mea)
plt.plot(f_ghz, rel_err)
plt.ylim((-20,10))

plt.tight_layout()
plt.show()


The agreement between measurements and the model seems very reasonnable. Relative error of $$\epsilon_{r,eff}$$ stay within ±1% outside very low frequencies and relative error of $$A$$ is kept between ±10% on most of the range. Considering the shape of $$A$$, it is not possible to do much better with this model.

### Connector effect characterization¶

The NISTMultilineTRL calibration coeffients contain information about the connector characteristics, which is corrected by the calibration.

These coefficients can be used to fit a connector model based on a transmission line section.

#### Delay and attenuation¶

In [14]:

# extract connector characteristic
# note: a half thru line is embedded in connector network
coefs = cal.coefs
r = skrf.mf.sqrt_phase_unwrap(coefs['forward reflection tracking'])
s1 = np.array([[coefs['forward directivity'],r],
[r, coefs['forward source match']]]).transpose()

conn_and_half_thru = line100mm.copy()
conn_and_half_thru.s = s1


A half thru line is embedded in the extracted network, we remove it with dembedding.

In [15]:

# half thru model, with already determined dielectric permittivity and tand.
half = m.line(50e-3, 'm', embed=True, z0=m.Z0_f)

conn = conn_and_half_thru ** half.inv
conn.name = 'Connector'


Estimate connector delay with linear regression on unwrapped phase.

In [16]:

# connector delay estimation by linear regression on the unwrapped phase
xlim = 9000 # used to avoid phase jump if any
phi_conn = (np.angle(conn.s[:xlim,1,0]))
z = np.polyfit(f[:xlim], phi_conn, 1)
p = np.poly1d(z)
delay_conn = -z[0]/(2*np.pi)
print('Connector delay: {:.1f} ps'.format(delay_conn * 1e12))

Connector delay: 42.5 ps


Build connector model and compare it against calibration extracted data.

In [17]:

mc = skrf.media.DefinedAEpTandZ0(m.frequency, ep_r=1, tanD=0.02, Z0=50,
f_low=1e3, f_high=1e18, f_ep=f_epr_tand, model='frequencyinvariant')

Z0_conn = 50.0 # the actual connector characteristic impedance will be tuned later
left = mc.line(delay_conn, 's', z0=Z0_conn)
check = mc.thru() ** left ** mc.thru()

plt.figure()
plt.subplot(2,1,1)
conn.plot_s_deg(1, 0, label='measured')
check.plot_s_deg(1, 0, label='model')
plt.legend()

plt.subplot(2,1,2)
conn.plot_s_db(1, 0, label='Measured')
check.plot_s_db(1, 0, label='Model')
plt.xlabel('Frequency (GHz)')
plt.ylabel('Insertion Loss (dB)')
plt.legend()

plt.tight_layout()
plt.show()

/home/docs/checkouts/readthedocs.org/user_builds/scikit-rf/conda/stable/lib/python3.5/site-packages/scikit_rf-0.14.9-py3.5.egg/skrf/mathFunctions.py:163: RuntimeWarning: divide by zero encountered in log10
out = 20 * npy.log10(input)


Comparison of connector model characteristics against calibration results shows a reasonnable agreement. Calibration results exhibit some glitches that do not correspond to the awaited physical behavior and may be related to singularities in the computations. It could possibly be enhanced by feeding more distinct lines to the algorithm, but these are not manufactured yet.

The results look good enough to be used.

#### Characteristic impedance¶

We now have estimated connector delay and attenuation, but what about the characterisic impedance? This value is required to properly parametrize the transmission line section model.

Optimization is used to find the characteristic impedance that minimize the difference between modelled and measured return loss.

In [18]:

s11_ref = conn.s[:,0,0]
x0 = [Z0_conn]

# function to be minimised
def model2(x, mc, delay_conn, s11_ref):
Z0_mod = x[0]
conn_mod = mc.line(delay_conn, 's', z0=Z0_mod)
check = mc.thru() ** conn_mod ** mc.thru()
s11_mod = check.s[:,0,0]

return sum(absolute(s11_ref-s11_mod))

# run optimizer
res = minimize(model2, x0, args=(mc, delay_conn, s11_ref),
bounds=[(50.5, 100)])

# get the results and print the results
Z0_conn = res.x[0]
print('Z0_conn={:.1f} ohm.'.format(Z0_conn))

Z0_conn=53.1 ohm.


The modelled results are plotted against the calibration data, as a sanity check.

In [19]:

conn_mod = mc.line(delay_conn, 's', z0=Z0_conn)
check = mc.thru() ** conn_mod ** mc.thru()

plt.figure;
plt.subplot(2,1,1)
check.plot_s_db(0,0)
conn.plot_s_db(0,0)

plt.subplot(2,1,2)
plt.plot(check.f*1e-9, (check.s_db[:,0,0]-conn.s_db[:,0,0]))
plt.ylabel('Delta (dB)')
plt.xlabel('Frequency (GHz)')

plt.tight_layout()
plt.show()


The delta in dB is quite big in low and high frequencies, but we will see when comparing measurement and simulation that the time domain reflectometry results are very decent with this value. The connector effects inductive peaks are well rendered in the case of the microstripline.

## Simulation¶

### Frequency axis¶

In [20]:

freq = skrf.F(1,10e3,10000, 'mhz')


### Media sections with different geometries¶

In [21]:

# 50 ohm segment
MSL1 = skrf.media.MLine(frequency=freq, z0=50, w=W, h=H, t=T,
ep_r=ep_r, mu_r=1, rho=1.712e-8, tand=tanD, rough=0.15e-6,
f_low=1e3, f_high=1e12, f_epr_tand=f_epr_tand,
diel='djordjevicsvensson', disp='kirschningjansen')

# Capacitive segment
MSL2 = skrf.media.MLine(frequency=freq, z0=50, w=8.0e-3, h=H, t=T,
ep_r=ep_r, mu_r=1, rho=1.712e-8, tand=tanD, rough=0.15e-6,
f_low=1e3, f_high=1e12, f_epr_tand=f_epr_tand,
diel='djordjevicsvensson', disp='kirschningjansen')

# Inductive segment
MSL3 = skrf.media.MLine(frequency=freq, z0=50, w=1.0e-3, h=H, t=T,
ep_r=ep_r, mu_r=1, rho=1.712e-8, tand=tanD, rough=0.15e-6,
f_low=1e3, f_high=1e12, f_epr_tand=f_epr_tand,
diel='djordjevicsvensson', disp='kirschningjansen')

# Connector transmission line media with guessed loss
MCON = skrf.media.DefinedAEpTandZ0(m.frequency, ep_r=1, tanD=0.025, Z0=50,
f_low=1e3, f_high=1e18, f_ep=f_epr_tand, model='djordjevicsvensson')


### Simulated devices under test¶

In [22]:

# SMA connector
conn  = MCON.line(delay_conn, 's', z0=Z0_conn, embed=True)

# microstripline
thru_simu = conn ** MSL1.line(100e-3, 'm', embed=True, z0=MSL1.Z0_f) ** conn
thru_simu.name = 'P1-MSL_Thru_100-P2'

# stepped impedance microstripline
step_simu = conn \
** MSL1.line(50e-3, 'm', embed=True, z0=MSL1.Z0_f) \
** MSL2.line(20e-3, 'm', embed=True, z0=MSL2.Z0_f) \
** MSL3.line(20e-3, 'm', embed=True, z0=MSL3.Z0_f) \
** MSL1.line(50e-3, 'm', embed=True, z0=MSL1.Z0_f) \
** conn
step_simu.name = 'P1-MSL_Stepped_140-P2'

# write simulated data to .snp files
write_data = False
if write_data:
step_simu.write_touchstone(dir='tdr_measurement_vs_simulation/simulation/')
thru_simu.write_touchstone(dir='tdr_measurement_vs_simulation/simulation/')


## Notes¶

### Time resolution¶

After DC point extrapolation, and neglecting the effect of windowing, the time resolution of measurement is :raw-latex:\begin{equation*} Resolution = \frac{1}{f_{span}} \end{equation*}

where $$Resolution$$ is the resolution in s and $$f_{span}$$ is the frequency span in Hz.

In our case, with $$f_{span} = 10 \mskip3mu\mathrm{[GHz]}$$ :raw-latex:\begin{equation*} Resolution \approx \frac{1}{10^{10}} \approx 100 \quad \mathrm{[ps]} \end{equation*}

With an effective relative dielectric permittivity of about 3.5, this give following physical resolution: :raw-latex:\begin{equation*} Resolution_{meter} = \frac{Resolution \dot{} c_0}{\sqrt{\epsilon{}_r}} \approx \frac{10^{-10} \dot{} 3\dot{}10^8}{\sqrt{3.5} } \approx 16 \quad \mathrm{[mm]} \end{equation*}

We will use reflexion measurement, so the actual range is divided by two, because the signal goes back and forth. Thus, the approximative distance resolution on our device under test will be 8mm. The discontinuities on the stepped impedance microstripline and the connector effects shall be well visible.

### Measurement range limitation to avoid alias response¶

The measurement range should be set in such a way the true response of the device is seen without repetitions (alias).

\begin{equation*} Range = \frac{c_0}{\sqrt{\epsilon{}_r} \dot{} \Delta{}_f} \end{equation*}

Where $$Range$$ is the range in m, $$c_0$$ the speed of light in m/s, $$\Delta{}_f$$ the frequency step in Hz and $$\epsilon{}_r$$ is the effective relative permittivity constant of the device under test.

With the given measurement setup and considering $$\Delta{}_f = 1 \mskip3mu\mathrm{[MHz]}$$, worst case $$\epsilon{}_r \approx 5.0$$, the range is :raw-latex:\begin{equation*} Range \approx \frac{3\dot{}10^8}{\sqrt{5} \dot{} 10^6} \approx 134 \quad \mathrm{[m]} \end{equation*} Relative dielectric permittivity for FR-4 is approximatively 5.0, the effective $$\epsilon{}_r$$ of microstripline geometry will be smaller because of the air on top, which give a safety margin, the range being underevaluated.

We will use reflexion measurement, so the actual range is divided by two, because the signal goes back and forth. However, our longest device being 200mm long, we can just quietly forget aliasing.

In [23]: