{"id":148,"date":"2025-05-24T06:59:16","date_gmt":"2025-05-24T06:59:16","guid":{"rendered":"https:\/\/synchrotron-light.net\/wp\/?page_id=148"},"modified":"2025-07-20T05:12:57","modified_gmt":"2025-07-20T05:12:57","slug":"chapter-4-radiation-emitted-by-charges-in-bending-magnets","status":"publish","type":"page","link":"https:\/\/synchrotron-light.net\/wp\/?page_id=148","title":{"rendered":"Chapter 4: Radiation emitted by charges in bending magnets"},"content":{"rendered":"<p>Below is a set of python codes associated with Chapter 4 of Daniele Pelliccia and David M. Paganin, &#8220;Synchrotron Light: A Physics Journey from Laboratory to Cosmos&#8221; (Oxford University Press, 2025).<\/p>\n<p>In order to run any of these python codes, you will need to include the following header file.<\/p>\n<pre class=\"theme:obsidian lang:python decode:true\">import numpy as np \r\nimport matplotlib.pyplot as plt\r\nfrom sys import stdout\r\nfrom mpl_toolkits.mplot3d import Axes3D\r\nfrom matplotlib import cm\r\nfrom matplotlib.colors import LightSource\r\nfrom scipy import special\r\n\r\nsynchrotron_style = {\r\n    'font.family': 'FreeSerif',\r\n    'font.size': 30,\r\n    'axes.linewidth': 2.0,\r\n    'axes.edgecolor': 'black',\r\n    'grid.color': '#5b5b5b',\r\n    'grid.linewidth': 1.2,\r\n}<\/pre>\n<h2 id=\"The-basics:-plotting-a-radial-field\">Synchrotron radiation for electrons and protons<\/h2>\n<p>See Fig. 4.2.<\/p>\n<pre class=\"theme:obsidian lang:python decode:true\">r_e = 2.817940323e-15 # m, classical electron radius\r\ne = 1.602176634e-19 # C, electron charge\r\nc = 2.99792458e8 # m\/s, speed of light \r\nm_e = 9.109383702e-31  # kg, electron mass\r\nEn_e = 0.510998950e-3 # GeV, electron rest energy. \r\n\r\nepsilon_0 = 8.854187813e-12  # F\/m, vacuum permittivity\r\nm_p = 1.672621924e-27 # kg, proton mass\r\nr_p = 1.0\/(4*np.pi*epsilon_0)*e**2\/(m_p*c**2) #m, the equivalent of r_e for protons\r\nEn_p = 0.93827208816 # GeV, proton rest energy\r\n\r\n# These constants are now in units of W T^{-2} GeV^{-2} \r\nC_e = (2\/3)*(r_e*e**2*c\/m_e)\/En_e**2\r\nC_p = (2\/3)*(r_p*e**2*c\/m_p)\/En_p**2<\/pre>\n<pre class=\"theme:obsidian lang:python decode:true\">B = 1.0\r\nEn = np.linspace(1, 10, 100)\r\nP_e = C_e*B**2*En**2\r\nP_p = C_p*B**2*En**2\r\n\r\nwith plt.style.context(('seaborn-v0_8-whitegrid')):\r\n    plt.rcParams.update(synchrotron_style)\r\n    fig = plt.figure(figsize=(9,8))\r\n    plt.semilogy(En, P_e, 'grey', lw = 4, label=\"Electron\")\r\n    plt.semilogy(En, P_p, 'k', lw = 4, label=\"Proton\")\r\n    \r\n    plt.legend(fontsize = 30, frameon=True, framealpha=0.9)\r\n    plt.xticks(fontsize = 30)\r\n    plt.yticks(fontsize = 30)\r\n    plt.xlabel('Accelerator energy (GeV)', fontsize=34, labelpad=10)\r\n    plt.ylabel(\"Emitted power (W)\", fontsize=34, labelpad=5)\r\n    plt.tight_layout()\r\n\r\nplt.show()<\/pre>\n<p><img loading=\"lazy\" decoding=\"async\" class=\"alignnone wp-image-149\" src=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.2.png\" alt=\"\" width=\"600\" height=\"533\" srcset=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.2.png 600w, https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.2-480x426.png 480w\" sizes=\"(min-width: 0px) and (max-width: 480px) 480px, (min-width: 481px) 600px, 100vw\" \/><\/p>\n<h2 style=\"text-align: left;\">Radio-frequency cavity<\/h2>\n<p>See Fig. 4.4.<\/p>\n<pre class=\"theme:obsidian lang:python decode:true\">omega = 2*np.pi\r\nphi = 1.2\r\nt = np.linspace(0, 1, 100)\r\nV = np.sin(omega*t+phi)\r\n\r\nwith plt.style.context(('seaborn-v0_8-whitegrid')):\r\n    plt.rcParams.update(synchrotron_style)\r\n    fig = plt.figure(figsize=(9,7))\r\n    \r\n    plt.plot(t, V, 'grey', lw = 4,zorder=2)\r\n    plt.scatter(t[26],V[26], s = 100, edgecolors='black',zorder=3 )\r\n    plt.scatter(t[23],V[23], marker = 'x', s = 100, c='black',linewidths=3, zorder=3 )\r\n    plt.scatter(t[29],V[29], marker = '+', s = 100, c='black',linewidths=3, zorder=3 )\r\n    \r\n    plt.plot([t[23], t[23]], [-1.1, V[23]], 'k--', zorder=1)\r\n    plt.plot([t[29], t[29]], [-1.1, V[29]], 'k--', zorder=1)\r\n    #plt.plot([t[26], t[26]], [-1.1, V[26]], 'k--', zorder=1)\r\n#     plt.legend(fontsize = 24, frameon=True)\r\n    plt.ylim(-1.1,1.1)\r\n    plt.xticks([0.0, 0.5, 1.], [\"0\", \"T\/2\", \"T\"], fontsize = 26)\r\n    #plt.xticks(fontsize = 24)\r\n    plt.yticks([-1, 0.0, 1.], [\"$-V_{0}$\", \"0\", \"$V_{0}$\"], fontsize = 26)\r\n    plt.xlabel('Time', fontsize=28, labelpad=10)\r\n    plt.ylabel(\"RF Voltage\", fontsize=28, labelpad=-10)\r\n\r\nplt.show()<\/pre>\n<p><img loading=\"lazy\" decoding=\"async\" class=\"alignnone wp-image-151\" src=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.4.png\" alt=\"\" width=\"600\" height=\"475\" srcset=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.4.png 600w, https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.4-480x380.png 480w\" sizes=\"(min-width: 0px) and (max-width: 480px) 480px, (min-width: 481px) 600px, 100vw\" \/><\/p>\n<h2 id=\"Simplified-plot-of-synchrotron-emission\">Simplified plot of synchrotron emission<\/h2>\n<p>The following code will generate one of the panels in Fig 4.8. Repeat the simulation by changing <span class=\"wp-katex-eq\" data-display=\"false\">\\beta<\/span>.<\/p>\n<pre class=\"theme:obsidian lang:python decode:true\">c = 2.99792458e8 # m\/s, speed of light \r\nbeta = 0.5\r\nrho = 10\r\nomega_0 = beta*c\/rho\r\ntheta = np.linspace(0, 2*np.pi, 100, endpoint=True)\r\nt = np.linspace(0.000001,0, 50)\r\n\r\nwith plt.style.context(('default')):\r\n    \r\n    plt.rcParams.update(synchrotron_style)\r\n    \r\n    fig = plt.figure(figsize=(6,6))\r\n    plt.plot(rho*np.cos(theta), rho*np.sin(theta), 'b', lw=2) \r\n    for i,j in enumerate(t):\r\n\r\n        cx = rho*np.cos(omega_0*j)\r\n        cy = rho*np.sin(omega_0*j)\r\n\r\n        a = cx + c*j*np.cos(theta)\r\n        b = cy + c*j*np.sin(theta)\r\n\r\n        plt.plot(b, a, 'k')\r\n\r\n    plt.xlabel(\"Distance (m)\", fontsize=28)\r\n    plt.ylabel(\"Distance (m)\", fontsize=28)\r\n    \r\n    plt.xticks(ticks = [-250, 0, 250],labels=(\"-250\", \"0\", \"250\"), fontsize = 20)\r\n    plt.yticks(ticks = [-250, 0, 250],labels=(\"-250\", \"0\", \"250\"), fontsize = 20)\r\n    \r\n    plt.show()<\/pre>\n<h1><img loading=\"lazy\" decoding=\"async\" class=\"alignnone size-full wp-image-153\" src=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.8.png\" alt=\"\" width=\"421\" height=\"393\" srcset=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.8.png 421w, https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.8-300x280.png 300w\" sizes=\"(max-width: 421px) 100vw, 421px\" \/><\/h1>\n<h2>Time-bandwidth uncertainty relation<\/h2>\n<p>See Fig. 4.11.<\/p>\n<pre class=\"theme:obsidian lang:python decode:true\">n_points = 1500\r\n\r\n# Time sequence\r\nt = np.linspace(-3, 3, n_points) # units of nanoseconds\r\ndt = t[2] - t[1]\r\n\r\nv_0 = 4 # Central frequency, 1\/ns\r\nomega_0 = 25 # 2*np.pi*4 # Central frequency\r\n# Frequency sequence\r\nv = (1.0\/ (n_points*dt))* (np.arange(n_points) - n_points\/\/2)\r\nomega = (2*np.pi\/ (n_points*dt))* (np.arange(n_points) - n_points\/\/2)\r\ndv = v[1]-v[0]\r\n\r\n# Pulse duration\r\ndelta_t1 = 0.6 #ns\r\ndelta_t2 = 0.08 #ns\r\ndelta_t3 = 0.04 #ns\r\n\r\n# Define pulse with Gaussian envelope\r\npulse1 =  np.cos(omega_0*t)  * np.exp(-t**2\/(2*delta_t1**2)) \r\npulse2 =  np.cos(omega_0*t)  * np.exp(-t**2\/(2*delta_t2**2)) \r\npulse3 =  np.cos(omega_0*t)  * np.exp(-t**2\/(2*delta_t3**2))\r\n\r\n# Fourier Transform of the pulse\r\nf_p1 = np.fft.fftshift(np.fft.fft(pulse1))\/np.sqrt(n_points)\r\nf_p2 = np.fft.fftshift(np.fft.fft(pulse2))\/np.sqrt(n_points)\r\nf_p3 = np.fft.fftshift(np.fft.fft(pulse3))\/np.sqrt(n_points)\r\n\r\ntick_size = 22\r\nfont_size = 30\r\nwith plt.style.context(('seaborn-v0_8-whitegrid')):\r\n    \r\n    plt.rcParams.update(synchrotron_style)\r\n    \r\n    fig, axs = plt.subplots(3,2, figsize=(16,16))\r\n    axs[0,0].plot(t, np.real(pulse1)**2, lw=4)\r\n    axs[0,0].tick_params(axis='x', labelsize=tick_size)\r\n    axs[0,0].tick_params(axis='y', labelsize=tick_size)\r\n    axs[0,0].set_xlabel('Time (ns)', fontsize=font_size)\r\n    axs[0,0].set_ylabel('Intensity', fontsize=font_size)\r\n    ###\r\n    axs[0,1].plot(omega, np.abs(f_p1)**2\/(np.abs(f_p1)**2).max(), lw=4)\r\n    axs[0,1].set_xlim([-100, 100])\r\n    axs[0,1].tick_params(axis='x', labelsize=tick_size)\r\n    axs[0,1].tick_params(axis='y', labelsize=tick_size)\r\n    axs[0,1].set_xlabel('Angular frequency (GHz)', fontsize=font_size)\r\n    axs[0,1].set_ylabel('Norm. power spectrum', fontsize=font_size)\r\n\r\n    axs[1,0].plot(t, np.real(pulse2)**2, lw=4)\r\n    axs[1,0].tick_params(axis='x', labelsize=tick_size)\r\n    axs[1,0].tick_params(axis='y', labelsize=tick_size)\r\n    axs[1,0].set_xlabel('Time (ns)', fontsize=font_size)\r\n    axs[1,0].set_ylabel('Intensity', fontsize=font_size)\r\n    ###\r\n    axs[1,1].plot(omega, np.abs(f_p2)**2\/(np.abs(f_p2)**2).max(), lw=4)\r\n    axs[1,1].set_xlim([-100, 100])\r\n    axs[1,1].tick_params(axis='x', labelsize=tick_size)\r\n    axs[1,1].tick_params(axis='y', labelsize=tick_size)\r\n    axs[1,1].set_xlabel('Angular frequency (GHz)', fontsize=font_size)\r\n    axs[1,1].set_ylabel('Norm. power spectrum', fontsize=font_size)\r\n\r\n    axs[2,0].plot(t, np.real(pulse3)**2, lw=4)\r\n    axs[2,0].tick_params(axis='x', labelsize=tick_size)\r\n    axs[2,0].tick_params(axis='y', labelsize=tick_size)\r\n    axs[2,0].set_xlabel('Time (ns)', fontsize=font_size)\r\n    axs[2,0].set_ylabel('Intensity', fontsize=font_size)\r\n    ###\r\n    axs[2,1].plot(omega, np.abs(f_p3)**2\/(np.abs(f_p3)**2).max(), lw=4)\r\n    axs[2,1].set_xlim([-100, 100])\r\n    axs[2,1].tick_params(axis='x', labelsize=tick_size)\r\n    axs[2,1].tick_params(axis='y', labelsize=tick_size)\r\n    axs[2,1].set_xlabel('Angular frequency (GHz)', fontsize=font_size)\r\n    axs[2,1].set_ylabel('Norm. power spectrum', fontsize=font_size)\r\n     \r\n    plt.tight_layout() #w_pad=2\r\n    \r\n    plt.show()<\/pre>\n<p><img loading=\"lazy\" decoding=\"async\" class=\"alignnone wp-image-154\" src=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.11-1014x1024.png\" alt=\"\" width=\"600\" height=\"606\" \/><\/p>\n<h2 id=\"An-infinity-of-stacked-carts\">Spectral plots<\/h2>\n<p>Radiated power as a function of the harmonic order in Fig. 4.14.<\/p>\n<pre class=\"theme:obsidian lang:python decode:true\">e = 1.602e-19  #C\r\nc = 3.0e8  #m\/s\r\nrho = 10 #m\r\nepsilon_0 = 8.854e-12  #F\/m\r\nP_0 = (1.0\/(8*np.pi**2*epsilon_0)) * (c*e**2 \/(rho**2)) \r\n\r\nn = np.linspace(1,int(5e12), int(1e7))\r\nwith plt.style.context(('seaborn-v0_8-whitegrid')):\r\n    \r\n    plt.rcParams.update(synchrotron_style)\r\n    \r\n    fig, ax = plt.subplots(2,1, figsize=(16,14))\r\n    \r\n    for th, theta in enumerate([0.0, 0.01]):\r\n       \r\n        for gamma in [100, 500, 3000]:    \r\n            beta = np.sqrt(1-1\/gamma**2)\r\n            power = P_0*n**2 *beta**4 *(special.jvp(n,n*beta*np.cos(theta))**2 + np.tan(theta)**2 * \\\r\n                         special.jv(n,n*beta*np.cos(theta))**2 \/ beta**2 )\r\n     \r\n            ax[th].loglog(n, power, lw = 4, label=\"$\\\\gamma=$\"+str(gamma))\r\n    \r\n        print(power.max()\/power[0])\r\n    \r\n    ax[0].legend(fontsize = 26, frameon=True, loc='lower left', framealpha=1)\r\n    ax[1].legend(fontsize = 26, frameon=True, loc='lower left', framealpha=1)\r\n    \r\n    ax[0].tick_params(axis='x', labelsize=24)\r\n    ax[0].tick_params(axis='y', labelsize=24)\r\n    ax[1].tick_params(axis='x', labelsize=24)\r\n    ax[1].tick_params(axis='y', labelsize=24)\r\n    \r\n    plt.xlabel('Harmonic order $n$', fontsize=28, labelpad=10)\r\n    \r\n    ax[0].set_ylabel(\"$\\\\frac{\\mathrm{d}P_{n}}{\\mathrm{d}\\Omega}$\", rotation=0, fontsize=32, labelpad=25)\r\n    ax[1].set_ylabel(\"$\\\\frac{\\mathrm{d}P_{n}}{\\mathrm{d}\\Omega}$\", rotation=0, fontsize=32, labelpad=25)\r\n    \r\n    ax[0].set_title(\"$\\\\theta=0$ mrad\", fontsize = 26, loc='left')\r\n    ax[1].set_title(\"$\\\\theta=10$ mrad\", fontsize = 26, loc='left')\r\n    \r\n    ax[0].set_ylim([1e-32, 1e-15])\r\n    ax[1].set_ylim([1e-32, 1e-15])\r\n    \r\n    plt.tight_layout() #h_pad=4\r\n    \r\nplt.show()<\/pre>\n<p><img loading=\"lazy\" decoding=\"async\" class=\"alignnone wp-image-156\" src=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.14-1024x890.png\" alt=\"\" width=\"600\" height=\"522\" \/><\/p>\n<h2>Plots of <span class=\"wp-katex-eq\" data-display=\"false\">J_n^{2}(n \\beta)<\/span> and <span class=\"wp-katex-eq\" data-display=\"false\">J_n^{&#039;2}(n \\beta)<\/span> vs <span class=\"wp-katex-eq\" data-display=\"false\">\\beta<\/span><\/h2>\n<p>See Fig. 4.15.<\/p>\n<pre class=\"theme:obsidian lang:python decode:true \">beta = np.linspace(0.97,0.9999999, 1000)\r\ntheta=0\r\nnn = 10000\r\nwith plt.style.context(('seaborn-v0_8-whitegrid')):\r\n    plt.rcParams.update(synchrotron_style)\r\n    fig = plt.figure(figsize=(9,7))\r\n    plt.semilogy(beta, special.jvp(nn,nn*beta*np.cos(theta))**2, 'k', lw = 4, label=\"$\\\\beta=$\"+str(beta))\r\n    plt.xlabel('$\\\\beta$', fontsize=28, labelpad=10)\r\n    plt.ylabel(\"$\\ J'_{n}^{2}(n \\\\beta \\\\cos \\\\theta)$\", rotation=90, fontsize=32, labelpad=25)\r\n    plt.xticks(ticks = [0.97, 0.98, 0.99, 1.0],fontsize = 30)\r\n    plt.yticks(fontsize = 30)<\/pre>\n<p><img loading=\"lazy\" decoding=\"async\" class=\"alignnone wp-image-158\" src=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.15.png\" alt=\"\" width=\"600\" height=\"435\" srcset=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.15.png 600w, https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.15-480x348.png 480w\" sizes=\"(min-width: 0px) and (max-width: 480px) 480px, (min-width: 481px) 600px, 100vw\" \/><\/p>\n<h2 id=\"2D-power-distribution\">2D map of the power distribution<\/h2>\n<p>See Fig. 4.16.<\/p>\n<pre class=\"theme:obsidian lang:python decode:true\">e = 1.602e-19  #C\r\nc = 3.0e8  #m\/s\r\nrho = 10 #m\r\nepsilon_0 = 8.854e-12  #F\/m\r\nP_0 = (1.0\/(8*np.pi**2*epsilon_0)) * (c*e**2 \/(rho**2)) \r\n\r\nnsteps = 350\r\nnmax = 10100000\r\nthsteps = 100\r\nn = np.arange(100,int((nmax-100)\/nsteps)*nsteps, int(nmax\/nsteps))\r\ntheta = np.linspace(-0.05,0.05, thsteps)\r\n\r\ngamma = 100\r\nbeta = np.sqrt(1-1\/gamma**2)#beta = 0.999\r\n\r\npower = np.zeros((nsteps, thsteps))\r\n\r\nfor i,j in enumerate(n):\r\n    for k, th in enumerate(theta):\r\n        \r\n        power[i,k] = P_0*j**2 *beta**4 *( \\\r\n                special.jvp(j,j*beta*np.cos(th))**2 + \\\r\n                np.tan(th)**2 * special.jv(j,j*beta*np.cos(th))**2 \/ beta**2 )\r\n\r\nx, y = np.meshgrid(n, 1000*theta)\r\n\r\nls = LightSource(180, 45)\r\n\r\nwith plt.style.context(('default')):\r\n    \r\n    plt.rcParams.update(synchrotron_style)\r\n    \r\n    fig, ax = plt.subplots(figsize=(12,8))\r\n\r\n    ctr = ax.contourf(x, y, power.T, 60, cmap=cm.viridis)\r\n    plt.xticks(fontsize = 32)\r\n    plt.yticks(fontsize = 32)\r\n    plt.xlabel('Harmonic order $n$ $(\\\\times 10^{-19})$', fontsize=26, labelpad=10)\r\n    plt.ylabel('$\\\\theta$ (mrad)', fontsize=26, labelpad=10)\r\n\r\n    plt.xticks(ticks = [2e6, 4e6, 6e6, 8e6, 10e6],labels=(\"2\", \"4\", \"6\", \"8\", \"10\"), fontsize = 32)\r\n\r\n    cbar = plt.colorbar(ctr)\r\n    cbar.set_ticks(ticks = [0, 0.4e-19, 0.8e-19, 1.2e-19],labels=(\"0.0\", \"0.4\", \"0.8\", \"1.2\"))\r\n    cbar.ax.tick_params(labelsize=32)\r\n\r\n    plt.show()<\/pre>\n<p><img loading=\"lazy\" decoding=\"async\" class=\"alignnone wp-image-157\" src=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.16-1024x732.png\" alt=\"\" width=\"600\" height=\"429\" \/><\/p>\n<h2 id=\"Longitudinal-Doppler-effect\">Asymptotic expansion of <span class=\"wp-katex-eq\" data-display=\"false\">J_n^{&#039;}(n \\beta)<\/span><\/h2>\n<p>See Fig. 4.17. Note: the labels are added separately in the published book version.<\/p>\n<pre class=\"theme:obsidian lang:python decode:true\">gamma = 800.0 \r\nbeta = np.sqrt(1.0 - 1.0\/gamma**2)\r\n\r\nnsteps = 5000\r\nnc = 1.5*gamma**3\r\nnmax = int(3*nc)\r\nn = np.arange(1,int((nmax-100)\/nsteps)*nsteps, int(nmax\/nsteps))\r\n\r\noriginal_function =  n**2 * special.jvp(n,n*beta)**2\r\nasymptotic_expansion = (n\/(2*np.pi*gamma)) * np.exp(-2*n\/(3*gamma**3))\r\n\r\nwith plt.style.context(('seaborn-v0_8-whitegrid')):\r\n    plt.rcParams.update(synchrotron_style)\r\n    fig, axs = plt.subplots(1, 2,figsize=(16,7))\r\n    fig.subplots_adjust(hspace=0)\r\n\r\n    # Note: the labels are added separately in the published book version\r\n    axs[0].plot(n\/nc, original_function, 'k', lw = 4, label=\"                    \")\r\n    axs[0].plot(n\/nc, asymptotic_expansion, 'gray', lw = 4, label=\"                \")\r\n    axs[0].set_ylabel('Comparison', fontsize=28, labelpad=10)\r\n    axs[0].tick_params(axis='y', labelsize=26)\r\n    axs[0].tick_params(axis='x', labelsize=26)\r\n    axs[0].set_xlabel('$n\/n_{c}$', fontsize=28, labelpad=10)\r\n        \r\n    axs[1].plot(n\/nc, (original_function-asymptotic_expansion)\/original_function, 'k', lw = 4)\r\n    axs[1].set_xlabel('$n\/n_{c}$', fontsize=28, labelpad=10)\r\n    axs[1].set_ylabel('Relative error', fontsize=28, labelpad=0)\r\n    axs[1].tick_params(axis='y', labelsize=26)\r\n    axs[1].tick_params(axis='x', labelsize=26)\r\n    \r\n    axs[0].legend(fontsize = 36, frameon=True, framealpha = 1, labelspacing=0.5, handlelength = 1, loc=(0.14,0.03))<\/pre>\n<p><img loading=\"lazy\" decoding=\"async\" class=\"alignnone size-full wp-image-160\" src=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.17.png\" alt=\"\" width=\"1010\" height=\"462\" srcset=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.17.png 1010w, https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.17-980x448.png 980w, https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.17-480x220.png 480w\" sizes=\"(min-width: 0px) and (max-width: 480px) 480px, (min-width: 481px) and (max-width: 980px) 980px, (min-width: 981px) 1010px, 100vw\" \/><\/p>\n<h2>Numerical computation of critical harmonic number<\/h2>\n<p>See Fig. 4.18.<\/p>\n<pre class=\"theme:obsidian lang:python decode:true\">e = 1.602e-19  #C\r\nc = 3.0e8  #m\/s\r\nrho = 10 #m\r\nepsilon_0 = 8.854e-12  #F\/m\r\nP_0 = (1.0\/(8*np.pi**2*epsilon_0)) * (c*e**2 \/(rho**2)) \r\n\r\nnsteps = 500\r\nnmax = 10000000\r\nthsteps = 501\r\n\r\ntheta = np.linspace(-0.05,0.05, thsteps)\r\ndth = theta[1]-theta[0]\r\n\r\nwith plt.style.context(('seaborn-whitegrid')):\r\n    \r\n    plt.rcParams.update(synchrotron_style)\r\n    \r\n    fig = plt.figure(figsize=(12,8))\r\n    for gamma in [10,20, 50, 500,2000]:\r\n        \r\n        nc = 1.5*gamma**3\r\n        \r\n        n = np.arange(1,5*nc, int(5*nc\/nsteps))\r\n        if n.shape[0]&gt;nsteps:\r\n            n = n[:nsteps]\r\n        \r\n        beta = np.sqrt(1-1\/gamma**2)#beta = 0.999\r\n\r\n        power2d = np.zeros((nsteps, thsteps))\r\n\r\n        for i,j in enumerate(n):\r\n            for k, th in enumerate(theta):\r\n\r\n                power2d[i,k] = (27*gamma**2\/(8*nc**2))*j**2 *( \\\r\n                        special.jvp(j,j*beta*np.cos(th))**2 + \\\r\n                        np.tan(th)**2 * special.jv(j,j*beta*np.cos(th))**2 \/ beta**2 )*np.cos(th)\r\n\r\n        power1d = dth*np.sum(power2d, axis=1)     \r\n\r\n        cumul = np.zeros_like(n).astype('float32')\r\n        for i,j in enumerate(n): \r\n            cumul[i]  = np.sum(power1d[:i])\/np.sum(power1d)  \r\n\r\n        plt.plot(n\/nc, cumul,  lw = 4,label=\"$\\\\gamma=$\"+str(gamma))\r\n    \r\n    \r\n    plt.yticks(ticks = [0, 0.25, 0.5, 0.75, 1], fontsize = 46)\r\n    plt.xticks(fontsize = 46)\r\n    plt.xlabel('$n\/n_{c}$', fontsize=32, labelpad=10)\r\n    plt.legend(fontsize = 36, frameon=True, loc='lower right', framealpha = 1, handlelength = 1)  \r\n    \r\n    plt.show()<\/pre>\n<p><img loading=\"lazy\" decoding=\"async\" class=\"alignnone wp-image-161\" src=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.18.png\" alt=\"\" width=\"600\" height=\"419\" srcset=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.18.png 600w, https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.18-480x335.png 480w\" sizes=\"(min-width: 0px) and (max-width: 480px) 480px, (min-width: 481px) 600px, 100vw\" \/><\/p>\n<h2>Modified Bessel functions of fractional order<\/h2>\n<p>See Fig. 4.19.<\/p>\n<pre class=\"theme:obsidian lang:python decode:true\">x = np.logspace(-5, 1, 10000)\r\n\r\nK13 = special.kv(1\/3,x)\r\nK23 = special.kv(2\/3,x)\r\n\r\nwith plt.style.context(('seaborn-v0_8-whitegrid')):\r\n    \r\n    plt.rcParams.update(synchrotron_style)\r\n        \r\n    fig, axs = plt.subplots(1, 1, figsize=(12,10))\r\n\r\n    axs.loglog(x, K13, 'k', lw = 4, label=\"$K_{1\/3}(z) $\")\r\n    axs.plot(x, K23, 'gray', lw = 4, label=\"$K_{2\/3}(z)$\")\r\n    axs.set_ylabel('Modified Bessel functions', fontsize=36, labelpad=0)\r\n    axs.set_xlabel('$z$', fontsize=36, labelpad=10)\r\n    axs.legend(fontsize = 36, frameon=True, labelspacing=1.0, loc='lower left', framealpha = 1, handlelength = 1)\r\n    axs.tick_params(axis='x', labelsize=46, pad=10)\r\n    axs.tick_params(axis='y', labelsize=46)\r\n    axs.yaxis.get_offset_text().set_fontsize(20)\r\n    \r\nplt.show()<\/pre>\n<p><img loading=\"lazy\" decoding=\"async\" class=\"alignnone wp-image-162\" src=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.19-1024x847.png\" alt=\"\" width=\"600\" height=\"496\" \/><\/p>\n<h2 id=\"Comparison-between-exact-and-approximated-spectrum\">Comparison between exact and approximated spectrum<\/h2>\n<p>See Fig. 4.20.<\/p>\n<pre class=\"theme:obsidian lang:python decode:true \">e = 1.602e-19  #C\r\nc = 3.0e8  #m\/s\r\nrho = 10 #m\r\nepsilon_0 = 8.854e-12  #F\/m\r\n\r\ntheta = 0.0000\r\nnsteps = 5000\r\n\r\nwith plt.style.context(('seaborn-v0_8-whitegrid')):\r\n    plt.rcParams.update(synchrotron_style)\r\n    fig, axs = plt.subplots(2,2, figsize=(16,14))\r\n    fig.subplots_adjust(hspace=0.3)\r\n    \r\n    for ax, gamma in zip(axs.flat, [2,10,100,1000]):    \r\n\r\n        beta = np.sqrt(1.0-1.0\/gamma**2)\r\n        \r\n        nc = 1.5*gamma**3\r\n        nmax = int(5.*nc)\r\n        n = np.linspace(1, nmax, nsteps)\r\n\r\n        P_tot = (1.0\/(6*np.pi*epsilon_0)) * (c*e**2 \/(rho**2)) * beta**4 * gamma**4 \r\n\r\n        arg_j = (n*beta*np.cos(theta))#.astype('uint64')\r\n        exact = P_tot*3*(n\/nc)**2 * gamma**2 *(special.jvp(n,arg_j)**2 + np.tan(theta)**2 * \\\r\n                                 special.jv(n,arg_j)**2 \/ beta**2 )\/np.pi\r\n        arg_k = 0.5*n\/nc * (1+theta**2 * gamma**2)**1.5\r\n        approximated = P_tot*(1.0\/np.pi**2)*(n\/nc)**2 * \\\r\n                       ( (1+theta**2 * gamma**2)**2 *special.kv(2\/3,arg_k)**2 \/gamma**2 + \\\r\n                                       theta**2 *(1+theta**2 * gamma**2)* special.kv(1\/3,arg_k) )\/np.pi\r\n\r\n        ax.plot(n\/nc, exact, 'k', lw = 4, label=\"Exact\")\r\n        ax.plot(n\/nc, approximated, 'gray', lw = 4, label=\"Approximated\")\r\n        ax.set_ylabel('$\\\\frac{\\mathrm{d}P_{n}}{\\mathrm{d}\\Omega}$', rotation=0, fontsize=32, labelpad=25)\r\n        ax.set_xlabel('$n\/n_{c}$', fontsize=28, labelpad=1)\r\n        ax.legend(fontsize = 24, frameon=True, labelspacing=0.4, loc='upper right')\r\n        ax.set_xticks([0,1,2,3,4,5], labels=['0','1','2','3','4','5'])\r\n        ax.tick_params(axis='x', labelsize=24)\r\n        ax.tick_params(axis='y', labelsize=24)\r\n        ax.set_title('$\\\\gamma = $'+str(int(gamma)), fontsize=28, loc='right', pad=10)\r\n        ax.yaxis.get_offset_text().set_fontsize(24)\r\n<\/pre>\n<p><img loading=\"lazy\" decoding=\"async\" class=\"alignnone size-full wp-image-164\" src=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.20.png\" alt=\"\" width=\"985\" height=\"863\" srcset=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.20.png 985w, https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.20-980x859.png 980w, https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.20-480x421.png 480w\" sizes=\"(min-width: 0px) and (max-width: 480px) 480px, (min-width: 481px) and (max-width: 980px) 980px, (min-width: 981px) 985px, 100vw\" \/><\/p>\n<h2>Surface plot of universal function<\/h2>\n<p>See Fig. 4.21.<\/p>\n<pre class=\"theme:obsidian lang:python decode:true \">e = 1.602e-19  #C\r\nc = 3.0e8  #m\/s\r\nrho = 10 #m\r\nepsilon_0 = 8.854e-12  #F\/m\r\n\r\ngamma = 1000\r\nbeta = np.sqrt(1.0-1.0\/gamma**2)\r\n\r\nnc = int(1.5*gamma**3)\r\nnsteps = 350\r\nnmax = int(5*nc)\r\nn = np.linspace(1, nmax, nsteps)\r\n\r\nthsteps = 200\r\ntheta = np.linspace(-0.005,0.005, thsteps)\r\n\r\nP_tot = (1.0\/(6*np.pi*epsilon_0)) * (c*e**2 \/(rho**2)) * beta**4 * gamma**4 \r\npower2d = np.zeros((nsteps, thsteps))\r\nfor i,j in enumerate(n):\r\n    for k, th in enumerate(theta):\r\n        \r\n        argum = 0.5*j\/nc * (1+th**2 * gamma**2)**1.5\r\n        \r\n        power2d[i,k] = P_tot*(1.0\/np.pi**2)*(j\/nc)**2 * \\\r\n                       ( (1+th**2 * gamma**2)**2 *special.kv(2\/3,argum)**2 \/gamma**2 + \\\r\n                                       th**2 *(1+th**2 * gamma**2)* special.kv(1\/3,argum) )\/np.pi\r\n\r\nX, Y = np.meshgrid(theta*gamma, n\/nc)\r\nwith plt.style.context(('seaborn-v0_8-whitegrid')):\r\n    \r\n    plt.rcParams.update(synchrotron_style)\r\n    fig = plt.figure(figsize=(15,15))\r\n    \r\n    ax1 = fig.add_subplot(111, projection='3d')     \r\n\r\n    surf1 = ax1.plot_surface(X, Y, power2d, rstride=1, cstride=1, cmap=cm.viridis,\r\n                           linewidth=0, antialiased=False)\r\n    ax1.view_init(20, 60)\r\n    \r\n    z1 = np.around(np.max(power2d)\/2, decimals=2)\r\n    z2 = np.around(np.max(power2d), decimals=2)\r\n    ax1.set_zticks([0,1e-17,2e-17])\r\n    ax1.set_zticklabels([0,1,2], fontsize=34)\r\n    ax1.set_xticks([-4, -2, 0, 2, 4])\r\n    ax1.set_xticklabels([-4, -2, 0, 2, 4], fontsize=34)\r\n    ax1.set_yticks([0, 1, 2, 3, 4, 5])\r\n    ax1.set_yticklabels([0, 1, 2, 3, 4, 5], fontsize=34)\r\n    \r\n    ax1.set_xlabel('$\\\\gamma \\, \\\\theta$ (rad)', fontsize=36, labelpad=30)\r\n    ax1.set_ylabel('$n\/n_{c}$', fontsize=36, labelpad=30)\r\nplt.show()<\/pre>\n<p><img loading=\"lazy\" decoding=\"async\" class=\"alignnone wp-image-165\" src=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.21.png\" alt=\"\" width=\"600\" height=\"600\" srcset=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.21.png 600w, https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.21-480x480.png 480w\" sizes=\"(min-width: 0px) and (max-width: 480px) 480px, (min-width: 481px) 600px, 100vw\" \/><\/p>\n<h2>Taylor expansion of <span class=\"wp-katex-eq\" data-display=\"false\">g(\\tau)<\/span><\/h2>\n<p>Calculate Taylor expansion recursively and plot Fig. 4.22.<\/p>\n<pre class=\"theme:obsidian lang:python decode:true\">def cn_recursive(n):\r\n    # Base case: c_0 = 0\r\n    if n == 0:\r\n        return 0\r\n    # Recursive case: c_n = c_{n-1}\/2 + 1\/2^(n-1)\r\n    else:\r\n        return 0.5*cn_recursive(n-1) + 1\/2**(n-1)\r\n\r\ndef derivative_K(z, N):\r\n    if N == 0:\r\n        return z*special.kv(2\/3,z\/2)\r\n    else:\r\n        return cn_recursive(N)*special.kvp(2\/3,z\/2,n=N-1) + (z\/2**N)*special.kvp(2\/3,z\/2,n=N)\r\n\r\n# Calculate the Taylor series up to fourth order\r\nord = 4\r\na = np.zeros(ord)\r\ntayl_exp = np.zeros_like(x)\r\nfor i in range(ord):\r\n    test = derivative_K(x, i)\r\n    a[i] = test[30]\/np.math.factorial(i)\r\n    tayl_exp += a[i]*(x-x[30])**i<\/pre>\n<pre class=\"theme:obsidian lang:python decode:true \"># Plot the Taylor expansion alongside the original function\r\nwith plt.style.context(('seaborn-v0_8-whitegrid')):\r\n    \r\n    plt.rcParams.update(synchrotron_style)\r\n    \r\n    fig, ax = plt.subplots(figsize=(10,10))\r\n    ax.plot(x, x*special.kv(2\/3,x\/2), 'k', lw=4, label = \"Exact\")\r\n    ax.plot(x, tayl_exp, 'r', lw=4, label=\"Taylor\")\r\n    ax.set_xlim([0.2, 1.8])\r\n    ax.set_ylim([1.02, 1.23])\r\n    plt.legend(fontsize = 30, frameon=True, framealpha=1)\r\n    \r\n    plt.xticks(fontsize = 36)\r\n    plt.yticks(fontsize = 36)\r\n\r\n    plt.xlabel('$\\\\tau$', fontsize=36)\r\n    plt.ylabel(\"$g(\\\\tau)$\", fontsize=36)<\/pre>\n<p><img loading=\"lazy\" decoding=\"async\" class=\"alignnone wp-image-167\" src=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.22.png\" alt=\"\" width=\"600\" height=\"525\" srcset=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.22.png 600w, https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.22-480x420.png 480w\" sizes=\"(min-width: 0px) and (max-width: 480px) 480px, (min-width: 481px) 600px, 100vw\" \/><\/p>\n<h2 id=\"Pade'-approximant\">Pad\u00e9 approximant<\/h2>\n<p>The coefficients of the Pad\u00e9 approximant are derived from the Taylor series coefficients using the procedure described in the book by Baker (1975) that is cited on page 134. Note that the coefficients a[i] calculated in the previous section are needed for the calculations described in the computational exercise on page 111 and to plot Fig. 4.23.<\/p>\n<pre class=\"theme:obsidian lang:python decode:true \">A = a[0]\r\nD = (a[2]**2 - a[3]*a[1]) \/ (a[1]**2 -a[2]*a[0])\r\nC = -a[0]*D\/a[1] - a[2]\/a[1]\r\nB = a[1]+a[0]*C\r\n\r\npade = (A+B*(x-1)) \/ (1+C*(x-1)+D*(x-1)**2)\r\n\r\nwith plt.style.context(('seaborn-v0_8-whitegrid')):\r\n    plt.rcParams.update(synchrotron_style)\r\n    fig, axs = plt.subplots(1,2, figsize=(21,8))\r\n    axs[0].plot(x, x*special.kv(2\/3,x\/2), 'k', lw=3,label = \"Exact\")\r\n    axs[0].plot(x, tayl_exp, 'r', lw=3,label=\"Taylor\")\r\n    axs[0].plot(x, pade, 'b', lw=3,label=\"Pade'\")\r\n    axs[0].set_xlim([0.18, 1.8])\r\n    axs[0].set_ylim([1, 1.25])\r\n    \r\n    axs[0].tick_params(axis='x', labelsize=32)\r\n    axs[0].tick_params(axis='y', labelsize=32)\r\n    axs[0].set_xlabel('$\\\\tau$', fontsize=36)\r\n    axs[0].set_ylabel(\"$g(\\\\tau)$\", fontsize=36)  \r\n    axs[0].legend(fontsize = 34, frameon=True, framealpha=0.9, loc = 'upper right', labelspacing=0.3, handlelength=1)\r\n    \r\n    axs[1].plot(x, x*special.kv(2\/3,x\/2), 'k', lw=3,label = \"Exact\")\r\n    axs[1].plot(x, tayl_exp, 'r', lw=3,label=\"Taylor\")\r\n    axs[1].plot(x, pade, 'b', lw=3,label=\"Pade'\")\r\n    axs[1].set_xlim([0.18, 10])\r\n    axs[1].set_ylim([0.05, 1.25])\r\n    \r\n    axs[1].tick_params(axis='x', labelsize=32)\r\n    axs[1].tick_params(axis='y', labelsize=32)\r\n    axs[1].set_xlabel('$\\\\tau$', fontsize=36)\r\n    axs[1].set_ylabel(\"$g(\\\\tau)$\", fontsize=36)  \r\n    axs[1].legend(fontsize = 34, frameon=True, framealpha=0.9, labelspacing=0.3, handlelength=1)\r\nplt.show()<\/pre>\n<p><img loading=\"lazy\" decoding=\"async\" class=\"alignnone size-large wp-image-168\" src=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.23-1024x417.png\" alt=\"\" width=\"1024\" height=\"417\" srcset=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.23-980x399.png 980w, https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.23-480x195.png 480w\" sizes=\"(min-width: 0px) and (max-width: 480px) 480px, (min-width: 481px) and (max-width: 980px) 980px, (min-width: 981px) 1024px, 100vw\" \/><\/p>\n<h2>Synchrotron spectrum plotted as a function of energy<\/h2>\n<p>See Fig. 4.24.<\/p>\n<pre class=\"theme:obsidian lang:python decode:true\">e = 1.602e-19  #C\r\nc = 3.0e8  #m\/s\r\nrho = 10 #m\r\nepsilon_0 = 8.854e-12  #F\/m\r\nhbar = 6.582119569e-16 #eV s\r\n\r\ntheta = 0.0000\r\nnsteps = 25000\r\n\r\nwith plt.style.context(('seaborn-v0_8-whitegrid')):\r\n    \r\n    plt.rcParams.update(synchrotron_style)\r\n    fig, ax = plt.subplots(1,1, figsize=(14,8))\r\n    \r\n    for gamma in [1000,5000,10000]:    \r\n\r\n        beta = np.sqrt(1.0-1.0\/gamma**2)\r\n        omega_0 = beta*c\/rho\r\n        \r\n        nc = 1.5*gamma**3\r\n        \r\n        n = np.logspace(8, 14, nsteps)\r\n        energy = hbar*omega_0*n\/1000 # energy units in keV\r\n        P_tot = (1.0\/(6*np.pi*epsilon_0)) * (c*e**2 \/(rho**2)) * beta**4 * gamma**4 \r\n\r\n        arg_k = 0.5*n\/nc * (1+theta**2 * gamma**2)**1.5\r\n        approximated = P_tot*(1.0\/np.pi**2)*(n\/nc)**2 * \\\r\n                       ( (1+theta**2 * gamma**2)**2 *special.kv(2\/3,arg_k)**2 \/gamma**2 + \\\r\n                                       theta**2 *(1+theta**2 * gamma**2)* special.kv(1\/3,arg_k) )\/np.pi\r\n\r\n        ax.loglog(energy, approximated, lw = 4, label=\"$\\\\gamma=$\"+str(gamma)+\"00\")\r\n        ax.set_ylabel('$\\\\frac{\\mathrm{d}^{2} \\\\mathcal{P}}{\\mathrm{d}\\Omega\\mathrm{d}\\\\mathcal{E}}$', rotation=0, fontsize=32, labelpad=35)\r\n        ax.set_xlabel('$\\\\mathcal{E} \\, (keV)$', fontsize=28, labelpad=5)\r\n        ax.legend(fontsize = 30, frameon=True, labelspacing=0.5, loc='lower left', framealpha=1, handlelength=1)\r\n        ax.tick_params(axis='x', labelsize=30, pad=10)\r\n        ax.tick_params(axis='y', labelsize=30)\r\n        \r\n        ax.set_ylim(1e-22, 4e-15)\r\n        ax.yaxis.get_offset_text().set_fontsize(20)\r\nplt.show()<\/pre>\n<p><img loading=\"lazy\" decoding=\"async\" class=\"alignnone wp-image-166\" src=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.24.png\" alt=\"\" width=\"600\" height=\"337\" srcset=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.24.png 600w, https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.24-480x270.png 480w\" sizes=\"(min-width: 0px) and (max-width: 480px) 480px, (min-width: 481px) 600px, 100vw\" \/><\/p>\n<h2>Polarisation: basic concepts<\/h2>\n<p>Plots of the polarisation ellipses in Fig. 4.27.<\/p>\n<pre class=\"theme:obsidian lang:python decode:true \">k_vec = 1\r\nz = np.linspace(0,10,101)\r\n# Chnage phi_x and phi_y to generate the other cases\r\nphi_x = 0\r\nphi_y = 3*np.pi\/4#\r\n\r\nEx_mod = 1\r\nEy_mod = 1\r\n\r\nEx = Ex_mod*np.cos(k_vec*z + phi_x)\r\nEy = Ey_mod*np.cos(k_vec*z + phi_y)\r\n\r\nwith plt.style.context(('seaborn-v0_8-whitegrid')):\r\n    plt.rcParams.update(synchrotron_style)\r\n    fig = plt.figure(figsize=(7,7))\r\n    plt.plot(Ex,Ey, 'k', lw = 4)\r\n    \r\n    plt.xticks(ticks = [-1,-0.5,0,0.5,1],fontsize = 30)\r\n    plt.yticks(ticks = [-1,-0.5,0,0.5,1], fontsize = 30)\r\n    plt.xlabel('x', fontsize=30)\r\n    plt.ylabel(\"y\", fontsize=30)\r\n    plt.title(\"$\\\\varphi_{x} = 0$, $\\\\varphi_{y} = 3\\pi\/4$\", fontsize=24)<\/pre>\n<p><img loading=\"lazy\" decoding=\"async\" class=\"alignnone size-full wp-image-169\" src=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.27.png\" alt=\"\" width=\"495\" height=\"488\" srcset=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.27.png 495w, https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.27-480x473.png 480w\" sizes=\"(min-width: 0px) and (max-width: 480px) 480px, (min-width: 481px) 495px, 100vw\" \/><\/p>\n<h2>3D plots of polarisation components<\/h2>\n<p>See Fig. 4.28.<\/p>\n<pre class=\"theme:obsidian lang:python decode:true \"># These are the values for the components [phi_x, phi_y, Ex_mod, Ey_mod] to generate the different cases\r\n# Linear 0: [0,0,1,0]\r\n# Linear 45: [0,0,1,1]\r\n# Circular right: [0,-pi\/2,1,1]\r\n# Circular left: [0,pi\/2,1,1]\r\n\r\nk_vec = 2*np.pi\/8\r\nz = np.linspace(0,20,101)\r\nphi_x = 0\r\nphi_y = 0\r\n\r\nEx_mod = 1\r\nEy_mod = 1\r\n\r\nEx = Ex_mod*np.cos(k_vec*z + phi_x)\r\nEy = Ey_mod*np.cos(k_vec*z + phi_y)\r\nL = np.sqrt(Ex**2 + Ey**2)\r\n\r\nwith plt.style.context(('seaborn-v0_8-whitegrid')):\r\n    fig = plt.figure(figsize=(15,10))\r\n    plt.rcParams.update(synchrotron_style)\r\n    ax = plt.subplot(projection='3d')\r\n    ax.set_proj_type('persp')\r\n\r\n    for i, j in enumerate(z):\r\n        ax.quiver(0, j, 0,  Ex[i], 0, Ey[i], linewidths=2, arrow_length_ratio=0.05)\r\n    \r\n    # Plot the axis\r\n    ax.quiver(0,0,0, 0,22,0, arrow_length_ratio=0.01, linewidths=4, colors='black')\r\n\r\n    ax.view_init(elev=10, azim=25)\r\n\r\n\r\n    ax.set_ylim3d(0, 22)\r\n    ax.set_xlim3d(-1, 1)\r\n    ax.set_zlim3d(-1, 1)\r\n\r\n    ax.set_zticks([-1,0,1])\r\n    ax.set_zticklabels([-1,0,1], fontsize=26)\r\n    plt.xticks([-1, 0, 1], fontsize = 26)\r\n    plt.yticks([0,5,10,15,20], fontsize = 26)\r\n  \r\n    plt.show()<\/pre>\n<p><img loading=\"lazy\" decoding=\"async\" class=\"alignnone wp-image-170\" src=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.28.png\" alt=\"\" width=\"600\" height=\"431\" srcset=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.28.png 600w, https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.28-480x345.png 480w\" sizes=\"(min-width: 0px) and (max-width: 480px) 480px, (min-width: 481px) 600px, 100vw\" \/><\/p>\n<h2 id=\"Polarisation-of-synchrotron-radiation\">Polarisation ellipses of synchrotron radiation<\/h2>\n<p>See Fig. 4.30.<\/p>\n<pre class=\"theme:obsidian lang:python decode:true \">gamma = 1000\r\nnc = int(1.5*gamma**3)\r\n\r\nphi = np.linspace(0,2*np.pi,360)\r\n\r\nbeta = 1.0 - 0.5\/gamma**2\r\n\r\nwith plt.style.context(('seaborn-v0_8-whitegrid')):\r\n\r\n    plt.rcParams.update(synchrotron_style)\r\n    fig, ax = plt.subplots(figsize=(10,10))\r\n    for th in [0, 0.00005, 0.0001, 0.0005]:\r\n        \r\n        E_par = np.cos(phi)\r\n        E_perp = (np.tan(th)\/beta)*special.jv(nc,nc*beta*np.cos(th))*np.sin(phi) \/ \\\r\n                  special.jvp(nc,nc*beta*np.cos(th))\r\n        \r\n        ax.plot(E_par,E_perp, lw = 4, label = '$\\\\theta=$'+str(1000*th)+\" mrad\")\r\n    \r\n    plt.xticks(ticks = [-1,-0.5,0,0.5,1],fontsize = 36)\r\n    plt.yticks(ticks = [-1,-0.5,0,0.5,1], fontsize = 36)\r\n    plt.xlabel('Parallel', fontsize=36)\r\n    plt.ylabel(\"Perpendicular\", fontsize=36)\r\n    plt.legend(fontsize = 30, frameon=True, ncol=2, loc='upper center', framealpha=1, handlelength = 1, handletextpad=0.5)\r\n    plt.title(\"$\\\\gamma = 1000$\", fontsize=30)\r\n    plt.show()<\/pre>\n<p><img loading=\"lazy\" decoding=\"async\" class=\"alignnone wp-image-171\" src=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.30.png\" alt=\"\" width=\"600\" height=\"587\" srcset=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.30.png 600w, https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.30-480x469.png 480w\" sizes=\"(min-width: 0px) and (max-width: 480px) 480px, (min-width: 481px) 600px, 100vw\" \/><\/p>\n<h2>Eccentricity<\/h2>\n<p>See Fig. 4.31.<\/p>\n<pre class=\"theme:obsidian lang:python decode:true \">e = 1.602e-19  #C\r\nc = 3.0e8  #m\/s\r\nrho = 10 #m\r\nepsilon_0 = 8.854e-12  #F\/m\r\n\r\ngamma = 1000\r\nbeta = np.sqrt(1.0-1.0\/gamma**2)\r\n\r\nnc = int(1.5*gamma**3)\r\nnsteps = 350\r\nnmax = int(5*nc)\r\nn = np.linspace(1, nmax, nsteps)\r\n\r\nthsteps = 200\r\ntheta = np.linspace(0,0.003, thsteps)\r\n\r\neccentricity = np.zeros((nsteps, thsteps))\r\nfor i,j in enumerate(n):\r\n    for k, th in enumerate(theta):\r\n        \r\n        argum = 0.5*j\/nc * (1+th**2 * gamma**2)**1.5\r\n        \r\n        eccentricity[i,k] = np.sqrt( 1.0 - (th**2 * gamma**2 \/ (1+ th**2 * gamma**2))*\\\r\n                       (special.kv(1\/3,argum)**2 \/special.kv(2\/3,argum)**2 ) )\r\n\r\nX, Y = np.meshgrid(theta*gamma, n\/nc)\r\n\r\nwith plt.style.context(('default')):\r\n\r\n    plt.rcParams.update(synchrotron_style)\r\n    fig, ax = plt.subplots(figsize=(12,10))\r\n\r\n    ctr = ax.contourf(X, Y, eccentricity, 14, cmap=cm.gist_earth)\r\n    ax.contour(X, Y, eccentricity, 14, colors='k', alpha=0.5)\r\n    plt.xticks(fontsize = 36)\r\n    plt.yticks(fontsize = 36)\r\n    plt.ylabel('$n\/n_{c}$', fontsize=32, labelpad=10)\r\n    plt.xlabel('$\\\\gamma \\\\vert \\\\theta \\\\vert$ (rad)', fontsize=32, labelpad=10)\r\n\r\n    cbar = plt.colorbar(ctr)\r\n    cbar.ax.tick_params(labelsize=36)\r\n\r\n    plt.show()<\/pre>\n<p><img loading=\"lazy\" decoding=\"async\" class=\"alignnone wp-image-172\" src=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.31.png\" alt=\"\" width=\"600\" height=\"547\" \/><\/p>\n<h2>Spectral intensity associated with the two polarisation components<\/h2>\n<p>See Fig. 4.34. First calculate the two components.<\/p>\n<pre class=\"theme:obsidian lang:python decode:true \">epsilon_0 = 8.854e-12  #F\/m\r\nrho = 10 #m\r\ngamma = 1000\r\nbeta = np.sqrt(1.0-1.0\/gamma**2)\r\nnc = int(1.5*gamma**3)\r\nnsteps = 350\r\nnmax = int(5*nc)\r\nthsteps = 200\r\nn = np.arange(100,int((nmax-100)\/nsteps)*nsteps, int(nmax\/nsteps))\r\ntheta = np.linspace(-0.005,0.005, thsteps)\r\n\r\nparallel = np.zeros((nsteps, thsteps))\r\nperpendicular = np.zeros((nsteps, thsteps))\r\nP_tot = (1.0\/(6*np.pi*epsilon_0)) * (c*e**2 \/(rho**2)) * beta**4 * gamma**4 \r\n# Neglect the factor P_tot\r\nfor i,j in enumerate(n):\r\n    for k, th in enumerate(theta):\r\n        \r\n        argum = 0.5*j\/nc * (1+th**2 * gamma**2)**1.5\r\n        parallel[i,k] = (1.0\/np.pi**2)*(j\/nc)**2 * (1+th**2 * gamma**2)**2 *special.kv(2\/3,argum)**2\/np.pi\r\n        perpendicular[i,k] = (1.0\/np.pi**2)*(j\/nc)**2 * th**2 * gamma**2 * (1+th**2 * gamma**2)* \\\r\n                              special.kv(1\/3,argum)\/np.pi\r\ntot_power = parallel + perpendicular<\/pre>\n<p>Then generate the plots.<\/p>\n<pre class=\"theme:obsidian lang:python decode:true\">X, Y = np.meshgrid(theta*gamma, n\/nc)\r\nwith plt.style.context(('seaborn-v0_8-whitegrid')):\r\n    plt.rcParams.update(synchrotron_style)\r\n    fig = plt.figure(figsize=(15,15))\r\n    \r\n    ax1 = fig.add_subplot(121, projection='3d')     \r\n    ax2 = fig.add_subplot(122, projection='3d') \r\n\r\n    surf1 = ax1.plot_surface(X, Y, parallel, rstride=1, cstride=1, cmap=cm.viridis,\r\n                           linewidth=0, antialiased=False)\r\n    ax1.view_init(20, 60)\r\n    \r\n    z1 = np.around(np.max(tot_power)\/2, decimals=2)\r\n    z2 = np.around(np.max(tot_power), decimals=2)\r\n\r\n    ax1.set_zticks([0,z1,z2])\r\n    ax1.set_zticklabels([0,z1,z2], fontsize=20)\r\n    ax1.set_zlim3d(0,z2)\r\n    ax1.zaxis.set_tick_params(pad=10)\r\n    ax1.set_xticks([-4, -2, 0, 2, 4])\r\n    ax1.set_xticklabels([-4, -2, 0, 2, 4], fontsize=22)\r\n    ax1.set_yticks([0, 1, 2, 3, 4, 5])\r\n    ax1.set_yticklabels([0, 1, 2, 3, 4, 5], fontsize=22)\r\n    \r\n    ax1.set_box_aspect((1, 1, 1.05))\r\n    \r\n    ax1.set_xlabel('$\\\\gamma \\, \\\\theta$ (rad)', fontsize=22, labelpad=15)\r\n    ax1.set_ylabel('$n\/n_{c}$', fontsize=22, labelpad=10)\r\n    ax1.set_title('Parallel', fontsize=28)\r\n    \r\n    \r\n    surf2 = ax2.plot_surface(X, Y, perpendicular, rstride=1, cstride=1, cmap=cm.viridis,\r\n                           linewidth=0, antialiased=False) \r\n    ax2.view_init(20, 60)\r\n    ax2.set_zticks([0,z1,z2])\r\n    ax2.set_zticklabels([0,z1,z2], fontsize=20)\r\n    ax2.zaxis.set_tick_params(pad=10)\r\n    ax2.set_zlim3d(0,z2)\r\n    ax2.set_xticks([-4, -2, 0, 2, 4])\r\n    ax2.set_xticklabels([-4, -2, 0, 2, 4], fontsize=22)\r\n    ax2.set_yticks([0, 1, 2, 3, 4, 5])\r\n    ax2.set_yticklabels([0, 1, 2, 3, 4, 5], fontsize=22)\r\n    ax2.set_box_aspect((1, 1, 1.05))\r\n    \r\n    ax2.set_xlabel('$\\\\gamma \\, \\\\theta$ (rad)', fontsize=22, labelpad=15)\r\n    ax2.set_ylabel('$n\/n_{c}$', fontsize=22, labelpad=10)\r\n    ax2.set_title('Perpendicular', fontsize=28)\r\nplt.savefig(\"polarisation_spectral_intensity.png\", dpi=150)       \r\nplt.show()<\/pre>\n<p><img loading=\"lazy\" decoding=\"async\" class=\"alignnone wp-image-173\" src=\"https:\/\/synchrotron-light.net\/wp\/wp-content\/uploads\/2025\/05\/shynchrotron-light-4.34-1024x545.png\" alt=\"\" width=\"600\" height=\"320\" \/><\/p>\n<p>&nbsp;<\/p>\n","protected":false},"excerpt":{"rendered":"<p>Below is a set of python codes associated with Chapter 4 of Daniele Pelliccia and David M. Paganin, &#8220;Synchrotron Light: A Physics Journey from Laboratory to Cosmos&#8221; (Oxford University Press, 2025). In order to run any of these python codes, you will need to include the following header file. import numpy as np import matplotlib.pyplot [&hellip;]<\/p>\n","protected":false},"author":1,"featured_media":0,"parent":79,"menu_order":0,"comment_status":"closed","ping_status":"closed","template":"","meta":{"_et_pb_use_builder":"","_et_pb_old_content":"","_et_gb_content_width":"","footnotes":""},"class_list":["post-148","page","type-page","status-publish","hentry"],"_links":{"self":[{"href":"https:\/\/synchrotron-light.net\/wp\/index.php?rest_route=\/wp\/v2\/pages\/148","targetHints":{"allow":["GET"]}}],"collection":[{"href":"https:\/\/synchrotron-light.net\/wp\/index.php?rest_route=\/wp\/v2\/pages"}],"about":[{"href":"https:\/\/synchrotron-light.net\/wp\/index.php?rest_route=\/wp\/v2\/types\/page"}],"author":[{"embeddable":true,"href":"https:\/\/synchrotron-light.net\/wp\/index.php?rest_route=\/wp\/v2\/users\/1"}],"replies":[{"embeddable":true,"href":"https:\/\/synchrotron-light.net\/wp\/index.php?rest_route=%2Fwp%2Fv2%2Fcomments&post=148"}],"version-history":[{"count":12,"href":"https:\/\/synchrotron-light.net\/wp\/index.php?rest_route=\/wp\/v2\/pages\/148\/revisions"}],"predecessor-version":[{"id":802,"href":"https:\/\/synchrotron-light.net\/wp\/index.php?rest_route=\/wp\/v2\/pages\/148\/revisions\/802"}],"up":[{"embeddable":true,"href":"https:\/\/synchrotron-light.net\/wp\/index.php?rest_route=\/wp\/v2\/pages\/79"}],"wp:attachment":[{"href":"https:\/\/synchrotron-light.net\/wp\/index.php?rest_route=%2Fwp%2Fv2%2Fmedia&parent=148"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}