{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "2fcad5e3",
   "metadata": {},
   "source": [
    "# Finding Dielectric Constant and Loss from Resonance Fitting"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7f62e543",
   "metadata": {},
   "source": [
    "In this example, the Q-factors are fitted from S-parameters data measured to characterize PCB laminate materials. The data correspond to a 36, 72 and 144 mm long stripline resonator respectively. We will then deduce the Dielectric permittivity and loss constants."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "89cc436a",
   "metadata": {},
   "source": [
    "## Fitting Q-factors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "cec79edf",
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "\n",
    "import skrf as rf\n",
    "\n",
    "rf.stylely()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d9862027",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Loading the data.\n",
    "# Setting the frequency unit is optional, just a convenience for displaying results\n",
    "Reson36 = rf.Network('data/resonator_36mm.s2p', f_unit='GHz')\n",
    "Reson72 = rf.Network('data/resonator_72mm.s2p', f_unit='GHz')\n",
    "Reson144 = rf.Network('data/resonator_144mm.s2p', f_unit='GHz')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9860ddc0",
   "metadata": {},
   "outputs": [],
   "source": [
    "Reson36.plot_s_db(m=1, n=0)\n",
    "Reson72.plot_s_db(m=1, n=0)\n",
    "Reson144.plot_s_db(m=1, n=0)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f77d1fbf",
   "metadata": {},
   "source": [
    "To find the Q-factors of these resonators, we create a `Qfactor` object from the transmission S-parameters we want to fit.\n",
    "\n",
    "Note that is there are many resonances close to each other, it is recommended to pass the frequency sub-Networks for the range of frequency of interest, to help the convergence of the fitting algorithm. Here, we are interested in the resonances occurring around 2 and 4 GHz, hence:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3e5f0990",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Around 2 GHz\n",
    "QF_36_2GHz = rf.Qfactor(Reson36['1.75-2.25GHz'].s21, res_type='transmission')\n",
    "QF_72_2GHz = rf.Qfactor(Reson72['1.75-2.25GHz'].s21, res_type='transmission')\n",
    "QF_144_2GHz = rf.Qfactor(Reson144['1.75-2.25GHz'].s21, res_type='transmission')\n",
    "# Around 4 GHz\n",
    "QF_36_4GHz = rf.Qfactor(Reson36['3.75-4.25GHz'].s21, res_type='transmission')\n",
    "QF_72_4GHz = rf.Qfactor(Reson72['3.75-4.25GHz'].s21, res_type='transmission')\n",
    "QF_144_4GHz = rf.Qfactor(Reson144['3.75-4.25GHz'].s21, res_type='transmission')\n",
    "\n",
    "all_QFs = [QF_36_2GHz, QF_72_2GHz, QF_144_2GHz, QF_36_4GHz, QF_72_4GHz, QF_144_4GHz]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "557dc5f1",
   "metadata": {},
   "outputs": [],
   "source": [
    "print(*all_QFs, sep='\\n')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2cc57ef0",
   "metadata": {},
   "source": [
    "Then we fit the Q-factors:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "aa771f3d",
   "metadata": {},
   "outputs": [],
   "source": [
    "[Q.fit() for Q in all_QFs]\n",
    "print(*all_QFs, sep='\\n')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2a27c7dd",
   "metadata": {},
   "source": [
    "We can create the Networks corresponding to the fitted results to benchmark the model against the measurements:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9e779273",
   "metadata": {},
   "outputs": [],
   "source": [
    "new_freq = rf.Frequency(1, 5, unit='GHz', npoints=1001)\n",
    "fitted_ntwk_36_2GHz = QF_36_2GHz.fitted_network(frequency=new_freq)\n",
    "fitted_ntwk_72_2GHz = QF_72_2GHz.fitted_network(frequency=new_freq)\n",
    "fitted_ntwk_144_2GHz = QF_144_2GHz.fitted_network(frequency=new_freq)\n",
    "\n",
    "fitted_ntwk_36_4GHz = QF_36_4GHz.fitted_network(frequency=new_freq)\n",
    "fitted_ntwk_72_4GHz = QF_72_4GHz.fitted_network(frequency=new_freq)\n",
    "fitted_ntwk_144_4GHz = QF_144_4GHz.fitted_network(frequency=new_freq)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c504a66b",
   "metadata": {},
   "outputs": [],
   "source": [
    "Reson36.plot_s_db(m=1, n=0, color='C0')\n",
    "fitted_ntwk_36_2GHz.plot_s_db(label='Fitted Model ~ 2GHz', lw=2, color='C0', ls='--')\n",
    "fitted_ntwk_36_4GHz.plot_s_db(label='Fitted Model ~ 4GHz', lw=2, color='C0', ls=':')\n",
    "Reson72.plot_s_db(m=1, n=0, color='C1')\n",
    "fitted_ntwk_72_2GHz.plot_s_db(label='Fitted Model ~ 2GHz', lw=2, color='C1', ls='--')\n",
    "fitted_ntwk_72_4GHz.plot_s_db(label='Fitted Model ~ 4GHz', lw=2, color='C1', ls=':')\n",
    "Reson144.plot_s_db(m=1, n=0, color='C2')\n",
    "fitted_ntwk_144_2GHz.plot_s_db(label='Fitted Model ~ 2GHz', lw=2, color='C2', ls='--')\n",
    "fitted_ntwk_144_4GHz.plot_s_db(label='Fitted Model ~ 4GHz', lw=2, color='C2', ls=':')\n",
    "plt.gca().legend(ncol=3)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "fed0d05b",
   "metadata": {},
   "source": [
    "Another way to represent the results is to use polar planes:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e4295d4e",
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(subplot_kw={'projection': 'polar'})\n",
    "Reson36.plot_s_polar(m=1, n=0, ax=ax, ls='', marker='x', ms=5)\n",
    "fitted_ntwk_36_2GHz.plot_s_polar(ax=ax, label=\"Fitted Model ~ 2GHz\", lw=2, ls='--', color='C3')\n",
    "fitted_ntwk_36_4GHz.plot_s_polar(ax=ax, label=\"Fitted Model ~ 4GHz\", lw=2, ls=':', color='C3')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "99f626cf",
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(subplot_kw={'projection': 'polar'})\n",
    "Reson72.plot_s_polar(m=1, n=0, ax=ax, ls='', marker='x', ms=5, color='C1')\n",
    "fitted_ntwk_72_2GHz.plot_s_polar(ax=ax, label=\"Fitted Model ~ 2GHz\", lw=2, ls='--', color='C3')\n",
    "fitted_ntwk_72_4GHz.plot_s_polar(ax=ax, label=\"Fitted Model ~ 4GHz\", lw=2, ls=':', color='C3')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "14f65cf4",
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(subplot_kw={'projection': 'polar'})\n",
    "Reson144.plot_s_polar(m=1, n=0, ax=ax, ls='', marker='x', ms=5, color='C2')\n",
    "fitted_ntwk_144_2GHz.plot_s_polar(ax=ax, label=\"Fitted Model ~ 2GHz\", lw=2, ls='--', color='C3')\n",
    "fitted_ntwk_144_4GHz.plot_s_polar(ax=ax, label=\"Fitted Model ~ 4GHz\", lw=2, ls=':', color='C3')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9f653dfd",
   "metadata": {},
   "source": [
    "## Increasing Accuracy\n",
    "If you want to get the best accuracy in the calculation of the Q-factors, it’s recommended to repeat the fitting a second time using S-parameters data having the resonance peak at the centre frequency and the span reduced to $\\pm f_L/Q_L$ (cf. [1] page 21). This is justified in two ways:\n",
    "\n",
    "* $Q_L$ is defined at the resonant frequency $f_L$\n",
    "* Limiting the span reduces the effects of other modes, mismatches, and differences between the actual resonator and the lumped-component model for the resonator. It also adds more points near the peak.\n",
    "\n",
    "The advantage is greatest if the resonance has poor shape, or if there are other modes nearby. Thus, a recommended workflow for optimal accuracy would be to:\n",
    "\n",
    "1. Set sweep parameters so the resonance is on the VNA screen\n",
    "2. Get Measurements\n",
    "3. Fit Q param to deduce $f_L$\n",
    "4. Set the Centre Freq to $f_L$ \n",
    "5. Set Span to $2f_L/Q_L$\n",
    "6. Get Measurements\n",
    "7. Fit Q param\n",
    "8. Calculate Unloaded Q\n",
    "\n",
    "So, for example:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2dc68e53",
   "metadata": {},
   "outputs": [],
   "source": [
    "# extracting the subfrequency range of interest from the data\n",
    "f_L = QF_72_2GHz.f_L\n",
    "Q_L = QF_72_2GHz.Q_L\n",
    "span = f_L/Q_L\n",
    "ntwk = Reson72[f'{f_L - span/2}-{f_L + span/2}Hz']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ce753704",
   "metadata": {},
   "outputs": [],
   "source": [
    "QF_72_2GHz_2 = rf.Qfactor(ntwk.s21, res_type='transmission')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "42c028f8",
   "metadata": {},
   "outputs": [],
   "source": [
    "# After the second fit, the weighting ratio is close to 0.2 at the minimum and maximum frequencies\n",
    "# provided that the loop_plan contains ‘w’. (cf. [1] Fig. 6(b) and equation (28))\n",
    "QF_72_2GHz_2.fit()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "53f9cc44",
   "metadata": {},
   "outputs": [],
   "source": [
    "QF_72_2GHz_2.Q_L  # higher accury fit"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "365f559c",
   "metadata": {},
   "source": [
    "\n",
    "\n",
    " "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ec60fb62",
   "metadata": {},
   "source": [
    "## Deducing Permittivity and Loss\n",
    "Equations for calculation Dk (aka $\\epsilon_r$) and Df (aka $\\tan \\delta$) from measured data for stripline resonators are given in [2]:\n",
    "\n",
    "$$\\epsilon_r = \\left[ \\frac{c n}{2 f_r (L + \\Delta L)} \\right]^2$$\n",
    "\n",
    "where:\n",
    "- $L$ is the length of resonator,\n",
    "- $n$ is the number of half wavelengths at resonance in resonator,\n",
    "- $f_r$ is the resonance frequency of resonator,\n",
    "- $c$ is the speed of light in vacuum,\n",
    "- $\\Delta L$ is the total effective increase in length of the resonant strip due to the fringing field at the ends of the resonant strip (neglected below for the example)\n",
    "\n",
    "and \n",
    "\n",
    "$$ \\tan \\delta = \\frac{1}{Q} - \\frac{1}{Q_c} $$\n",
    "\n",
    "where:\n",
    "- $Q_c$ is the Q-value associated with the copper loss: 250 at 2 GHz and 360 at 4 GHz according to IPC."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "50b26549",
   "metadata": {},
   "source": [
    "Hence, for the permittivity:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "396cab77",
   "metadata": {},
   "outputs": [],
   "source": [
    "def Dk(f, L, n):\n",
    "    \"Calculates Dk from resonator data. Eq.(1) of ref [1]\"\n",
    "    return (n*rf.c/(2*f*L))**2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4de8e99c",
   "metadata": {},
   "outputs": [],
   "source": [
    "Dk(QF_72_2GHz.f_L, 72e-3, 2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a3541406",
   "metadata": {},
   "outputs": [],
   "source": [
    "Dk(QF_72_2GHz.f_L, 144e-3, 4)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ab205dbe",
   "metadata": {},
   "source": [
    "And the loss tangent:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "62c94424",
   "metadata": {},
   "outputs": [],
   "source": [
    "Qc = 250  # from IPC ref [1]\n",
    "print(1/QF_36_2GHz.Q_L - 1/Qc)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "95609723",
   "metadata": {},
   "outputs": [],
   "source": [
    "Qc = 250  # from IPC ref [1]\n",
    "print(1/QF_144_2GHz.Q_L - 1/Qc)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "498648e8",
   "metadata": {},
   "outputs": [],
   "source": [
    "Qc = 360  # from IPC\n",
    "print(1/QF_36_4GHz.Q_L - 1/Qc)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e2c1e19a",
   "metadata": {},
   "outputs": [],
   "source": [
    "Qc = 360  # from IPC\n",
    "print(1/QF_144_4GHz.Q_L - 1/Qc)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9952f148",
   "metadata": {},
   "source": [
    "## References\n",
    "- [1] A. P. Gregory, \"Q-factor Measurement by using a Vector Network Analyser\", National Physical Laboratory Report MAT 58 (2021), https://eprintspublications.npl.co.uk/9304/\n",
    "- [2] IPC-TM-650 TEST METHODS MANUAL, https://www.ipc.org/sites/default/files/test_methods_docs/2-5_2-5-5-5.pdf"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "869021ad",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.11"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}