{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "\n# Simple dipole magnet\n\nBasic example of computation of the field emitted by an electron beam \ndeflected by a simple dipole magnet.\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Import modules required for simulations and for comparison with the \ntheoretical distribution of dipole radiation.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "import matplotlib.pyplot as plt\nimport numpy as np\nfrom scipy.special import kv\nfrom scipy.constants import e, epsilon_0, mu_0, c, pi, h\n\nimport pysrw as srw" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We take as example the parameters of a typical bending magnet of a \nthird-generation light source (ALBA-Spain)\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "energy = 3.0 # GeV\nrho = 7.047 # m\nlength = 1.384 # m\ngap = 36e-3 # m" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "which produce a particle deflection of\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "delta = 2 * np.arcsin(length / 2 / rho)\nprint(f\"Deflection angle: {delta*1e3:.2f} mrad\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The magnet is created with the default SRW model for dipoles \n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "dipole = srw.magnets.Dipole(energy=energy, bendingR=rho, \n coreL=length, edgeL=gap)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This dipole emits a broad band radiation with a critical wavelength \n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "wl_c = dipole.getLambdaCritical()\nprint(f\"Critical wavelength: {wl_c:.2f} nm\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The magnetic element must be included in a magnetic container.\nIn this simple example, the `dipole` is the only device and the limits for \nthe numerical integration as left as default.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "mag_container = srw.magnets.MagnetsContainer([dipole])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The last input for the simulation is the emitter, instance of \n:py:func:`~pysrw.emitters.ParticleBeam`\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "beam = srw.emitters.ParticleBeam(energy, xPos=0, yPos=0, zPos=0)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Before computing the wavefront, it is a good practice to obtain and plot \nthe particle trajectory and check that the simulation objects properly \ndefined\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "traj = srw.computeTrajectory(beam, mag_container)\nsrw.plotTrajectory(traj)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Create a 10 mm x 10 mm observation mesh, with a resolution of 100 um, \nplaced 1 m downstream of the source\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "observer = srw.wavefronts.Observer(centerCoord=[0, 0, 5],\n obsXextension=5e-3, \n obsYextension=5e-3,\n obsXres=200e-6,\n obsYres=20e-6)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can finally simulate the wavefront and derive the optical intensity,\nfor example at the critical wavelength. Note that, dealing with a dipole, \nwe use the \"auto-wigggler\" value for :py:func:`~pysrw.configuration.SR_APPROX`\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "wl = wl_c\nwfr = srw.computeSrWfrMultiProcess(4, particleBeam=beam, \n magnetsContainer=mag_container,\n observer=observer, wavelength=wl,\n srApprox=\"auto-wiggler\")\nintensity = wfr.getWfrI()\nsrw.plotI(intensity)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The distribution of SR emitted by a particle in a uniform field has an \nanalytical expression. Let's compare it to the simulation result.\nWe extract a vertical slice from the simulation radiation\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "cuts = srw.extractCrossCuts(intensity)\nyax = cuts[\"yax\"]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "and we derive the theoretical expression (Schwinger distributions)\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "r_p = observer.zPos\ntheta = np.arctan(yax / r_p)\ngamma = beam.nominalParticle.gamma\n\ngt_sqr = (gamma * theta)**2\nxi = wl_c/ (2 * wl) * (1 + gt_sqr)**(3/2)\n\ne_theo_h = (-np.sqrt(3)*e*gamma) / ((2*pi)**(3/2) * epsilon_0 * c *r_p) *\\\n (wl_c / (2 * wl)) * (1 + gt_sqr) * kv(2/3, xi)\n\ne_theo_v = (1j*np.sqrt(3)*e*gamma) / ((2 * pi)**(3/2) * epsilon_0 * c *r_p) *\\\n (wl_c / (2 * wl)) * (gamma * theta) * np.sqrt(1 + gt_sqr) * kv(1/3, xi)\n\nint_theo = (4 * pi) / (mu_0 * c * wl*1e-9 * h * e) \\\n * 1e-15 * (np.abs(e_theo_h)**2 + np.abs(e_theo_v)**2)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Simulated and theoretical distributions are finally plotted\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "fig, ax = plt.subplots()\nax.plot(theta*1e3, cuts[\"ydata\"])\nax.plot(theta*1e3, int_theo, label=\"theo\", linestyle=\"--\")\nax.set_xlabel(\"Polar angle [mrad]\")\nax.set_ylabel(\"Intensity [ph / (s mm^2 nm)]\")\nax.legend()\nplt.show()" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "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.9.9" } }, "nbformat": 4, "nbformat_minor": 0 }