{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Average radial distribution functions\n", "\n", "Here we calculate the average radial cumulative distribution functions between two groups of atoms.\n", "\n", "**Last executed:** Feb 06, 2020 with MDAnalysis 0.20.2-dev0\n", "\n", "**Last updated:** February 2020\n", "\n", "**Minimum version of MDAnalysis:** 0.17.0\n", "\n", "**Packages required:**\n", " \n", "* MDAnalysis (Michaud-Agrawal *et al.*, 2011, Gowers *et al.*, 2016)\n", "* MDAnalysisTests\n", "\n", "**Optional packages for visualisation:**\n", "\n", "* [matplotlib](https://matplotlib.org)\n" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import MDAnalysis as mda\n", "from MDAnalysis.tests.datafiles import TPR, XTC\n", "from MDAnalysis.analysis import rdf\n", "import matplotlib.pyplot as plt\n", "%matplotlib inline" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Loading files" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The test files we will be working with here feature adenylate kinase (AdK), a phosophotransferase enzyme. [[3]](#References)" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "u = mda.Universe(TPR, XTC)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Calculating the average radial distribution function for two groups of atoms" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "A radial distribution function $g_{ab}(r)$ describes the time-averaged density of particles in $b$ from the reference group $a$ at distance $r$. It is normalised so that it becomes 1 for large separations in a homogenous system.\n", "\n", "$$ g_{ab}(r) = (N_{a} N_{b})^{-1} \\sum_{i=1}^{N_a} \\sum_{j=1}^{N_b} \\langle \\delta(|\\mathbf{r}_i - \\mathbf{r}_j| - r) \\rangle$$\n", "\n", "The radial cumulative distribution function is \n", "\n", "$$G_{ab}(r) = \\int_0^r \\!\\!dr' 4\\pi r'^2 g_{ab}(r')$$.\n", "\n", "The average number of $b$ particles within radius $r$ at density $\\rho$ is:\n", "\n", "$$N_{ab}(r) = \\rho G_{ab}(r)$$\n", "\n", "The average number of particles can be used to compute coordination numbers, such as the number of neighbours in the first solvation shell.\n", "\n", "Below, I calculate the average RDF between each atom of residue 60 to each atom of water to look at the distribution of water over the trajectory. The RDF is limited to a spherical shell around each atom in residue 60 by `range`. Note that the range is defined around *each atom*, rather than the center-of-mass of the entire group.\n", "\n", "If you are after non-averaged radial distribution functions, have a look at the [site-specific RDF class](site_specific_rdf.ipynb). [The API docs for the InterRDF class are here.](https://docs.mdanalysis.org/stable/documentation_pages/analysis/rdf.html#MDAnalysis.analysis.rdf.InterRDF)" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "res60 = u.select_atoms('resid 60')\n", "water = u.select_atoms('resname SOL')\n", "\n", "irdf = rdf.InterRDF(res60, water,\n", " nbins=75, # default\n", " range=(0.0, 15.0), # distance in angstroms\n", " )\n", "irdf.run()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The distance bins are available at `irdf.bins` and the radial distribution function is at `irdf.rdf`. " ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([ 0.1, 0.3, 0.5, 0.7, 0.9, 1.1, 1.3, 1.5, 1.7, 1.9, 2.1,\n", " 2.3, 2.5, 2.7, 2.9, 3.1, 3.3, 3.5, 3.7, 3.9, 4.1, 4.3,\n", " 4.5, 4.7, 4.9, 5.1, 5.3, 5.5, 5.7, 5.9, 6.1, 6.3, 6.5,\n", " 6.7, 6.9, 7.1, 7.3, 7.5, 7.7, 7.9, 8.1, 8.3, 8.5, 8.7,\n", " 8.9, 9.1, 9.3, 9.5, 9.7, 9.9, 10.1, 10.3, 10.5, 10.7, 10.9,\n", " 11.1, 11.3, 11.5, 11.7, 11.9, 12.1, 12.3, 12.5, 12.7, 12.9, 13.1,\n", " 13.3, 13.5, 13.7, 13.9, 14.1, 14.3, 14.5, 14.7, 14.9])" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "irdf.bins" ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "scrolled": false }, "outputs": [ { "data": { "text/plain": [ "Text(0, 0.5, 'Radial distribution')" ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEGCAYAAABo25JHAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3deXxU9dX48c/JvidAICwhrAEElC2yue9YFaxa666tlS5atba12sXHatufbZ9a7aPV4lKsVXFXVBR3BEFllX3fEiAhELKQPZnz+2NuwiSZTCaYm0nIeb9eeWXuvd+5c4Zlznx3UVWMMcZ0XWGhDsAYY0xoWSIwxpguzhKBMcZ0cZYIjDGmi7NEYIwxXVxEqANordTUVB04cGCowzDGmE5l+fLlB1S1p79rnS4RDBw4kGXLloU6DGOM6VREZFdz16xpyBhjujhLBMYY08VZIjDGmC7OEoExxnRxlgiMMaaLs0RgjDFdnCUCY4zp4iwRGGNMB6eq/PGd9WzKLXHl/pYIjDGmg5u/Lo8nFu5g7Z4iV+5vicAYYzowj0d56MPNDEqNZ8bYvq68hiUCY4zpwN5fn8vG3BJuPWsoEeHufGRbIjDGmHb2/rpczvv7Z8z+fAcV1bXNlvPWBrYwODWe6WP6uRaPJQJjjGlH5VW1/M/cdewuKOPet9Zz8p8/5vEF2zhcWdOk7Hvr6moDmYSHiWsxWSIwxph29NSi7ewrqmD2907kxZmTOa5PEg+8u5GT//wxL3y1G49HAW9t4OEPtzCkZzwXjXGnb6BOp1uG2hhjOqv9JRU89uk2zhuVxqTBPQCYNLgHX2cX8sd5G7j7tTW8tiKH/3fJ8WzMLWFTXgkPXzHW1doAWCIwxph28/cPNlNZ4+Gu849rcH5M/xRenDmZl5fl8Md5Gzj/4YUkx0YytFcCF57gbm0ArGnIGGPaxcbcYl5cms11UwYyKDW+yXUR4fIT+/PRz0/jguP7cOBwFb84d5jrtQFwuUYgItOAh4Fw4ElVfaDR9b8DZziHcUAvVU1xMyZjjGlv3pnBG0iMieTWs4YGLJuaEM1DV4zj3umjSImLapf4XEsEIhIOPAqcA+QAS0Vkrqquryujqj/zKf9TYJxb8RhjTHspKK3iqx0FFJZVcaismpxDZSzccoDfXTgy6A/39koC4G6NYCKwVVW3A4jIHGAGsL6Z8lcC/+NiPMYY47qdB0q56okv2FtUUX8uKiKMk4b24NrJA0IYWfPcTAT9gGyf4xxgkr+CIjIAGAR83Mz1mcBMgIyMjLaN0hhj2sjW/Ye56okvqPEoz944kcE9E+gWF0lsZDgi7rf1H62OMmroCuAVVfU7xU5VZwGzALKysrQ9AzPGmGBsyi3h6ie/BOCFmyYzvHdiiCMKnpujhvYA/X2O051z/lwBvOBiLMYYE5TSyhruf3s9W/KCX/J5/d5irnziC8IE5szsXEkA3E0ES4FMERkkIlF4P+znNi4kIiOAbsASF2MxxpgWqSp3v7aGpxbt4IfPLqfUz7IPja3KLuTKJ74gJiKMl344haG9Etoh0rblWiJQ1RrgFmA+sAF4SVXXich9IjLdp+gVwBxVtSYfY0xIPfvFLuZ+vZeLxvRl58FSfvfm2oDlv9x+kGue/JLk2Ehe/OEUBvqZH9AZuNpHoKrzgHmNzt3T6PheN2MwxphgrNx9iPvfXs9ZI3rx8HfHMjg1noc/2sJJQ1K5dEJ6k/Kfbc5n5rPL6JcSy3M/mEzv5JgQRN02bGaxMabLKyit4ubnVpCWFMODl48lLEy49axMJg3qzu/eXMu2/MP1ZStranltRQ4/eGYZg1ITePGHUzp1EoCOM2rIGGPanaqSc6icX7++hgOHq3j1x1NJjosEIDxMePiKcZz/8Gfc8vxKfnDyID7amMeCTfmUVtUypn8Kz3zvxHad+OUWSwTGmC6ltLKG11fu4csdBSzdUUBusXfi1/+75HiOT09uULZ3cgx/u3wM35+9jJ+//DVpSdFMH9uPs0b04pRhqURHhIfiLbQ5SwTGmC7ji+0H+eUrX5NdUE5aUjQnDuzOpEHdmTy4B5lp/od8njkijed/MInEmEhG90vq0BPDjpYlAmPMMaOiupY/zdvAl9sLOHVYKtNG92Zc/25U1NTy53c38sySXQzsEcecmZOZNKh70B/qU4emuhx5aFkiMMYcE3YdLOXH/13B+n3FTBjQjdmLd/LEwh30TIwmKjyMPYXlfO+kgdx53ghio46NJp22YonAGNPpvbd2H798eTVhYcLTN2Rx5og0iiuq+WTjfuavyyWvuJIHLx9TvyuYacgSgTGmU3ty4Xb+8M4GxvRP4dGrxpHeLQ6ApJhIZoztx4yx/UIcYcdnicAY02kVlVfz0IdbOGN4T/51bRZRETY16mjYn5oxptN6/svdHK6s4RfnDbck8A3Yn5wxplOqrKnl35/v4JTMVEb1TW75CaZZlgiMMZ3Smyv3sr+kkh+eOiTUoXR6lgiMMZ2Ox6P867NtjOyTxElDbSTQN2WJwBjT6Xy0cT/b8kv54WmDj8mZvu3NRg0ZY1xTUFrF++tyuXRCOpHhLX/vrKn1MHvxTsqqaunfPZb0bnH07xZHWlJ0gw/8fy3YRr+UWC44vo+b4XcZlgiMMc2q9ShhwlF9666sqeUHzyxlxe5CdheUcee0EQHL19R6uP3FVby9el+Ta72TYjh9eE/OGNGL2Mhwlu06xP9cNJKIIJKLaZklAmOMX1U1Hi57fDFDeibw9++ObdVzVZXfvL6WFbsLGZ+RwmMLtjFlSA9Oyezpt3x1rYfb5qxk3ppc7j5/BNdPHciewnKyC8rYXVDGF9sP8s7qfcxZmg1ASlwk3z2xv997mdZzNRGIyDTgYSAceFJVH/BT5nLgXkCBr1X1KjdjMqYrW77rEGPSk4P6Jv3Yp9tYnVPE+r3F/PaC4+iREB306zy1aAevLM/htrMy+dFpQ5j+yCJ+9uLXvHvbKfRMbHifqhoPP31hBfPX5fHbC47jB6cMBmBIzwSG9PTu/3vdlIFU13pYvusQCzbnMyY9hbgo+x7bVlyrV4lIOPAocD4wErhSREY2KpMJ3A2cpKqjgNvdiseYrm51TiGXPraYxxdsa7HsptwSHvlkC1kDulHjUd5ctTfo1/l0037+NG8D54/uzW1nZRIbFc4jV42npKKaO15ahcdzZHvy7fmH+fF/lzN/XR73XDiyPgn4ExkexuTBPfjVtBFMG9076HhMy9xsYJsIbFXV7apaBcwBZjQqcxPwqKoeAlDV/S7GY0yX9u7aXACeXLSDw5U1zZarqfVw5ytfkxQTyazrsjghPZmXl+cEvHdVjYe1e4qY89VufvrCSob3TuJvl48hLMzbtzC8dyL3Th/Fwi0H+MfHW3hleQ6XP76EM/+2gE8353PfjFF8/+RBbfdmTau4WbfqB2T7HOcAkxqVGQYgIp/jbT66V1Xfa3wjEZkJzATIyMhwJVhjjmWqyntrc+nfPZbsgnKeXbKLH5/ufyLWU4t28HVOEY9cNY7u8VF8Z0I6v3tzHev2FjWZwfvRhjwe/GAzm/NKqK71ftPvlxLLE9dNaNJ0c8WJ/Vm09QAPfbgFgEGp8dw5bTiXjU+nV1Ln3vO3swt1I1sEkAmcDqQDn4nI8apa6FtIVWcBswCysrK08U2MMYFt2X+YHQdKuf/i0Xy4Po8nFm7n+qkDmnxYb88/zIMfbObckWn1QzMvGtOX+9/ewMvLchg1/UgiyC+p5GcvriI1IZobTx7MqL5JjO6XzIDucfU1AV8iwgOXHM+wXolMHtydia3YGMa4y82moT2Ab7d+unPOVw4wV1WrVXUHsBlvYjDGtKH5a3MRgfNGpnHrWZkUlFbx3Be7G5Qprazh5y9/TXREGH+4eHT9h3RKXBTnjErjzVV7qKrx1Jf/wzvrqaj2MOu6LO46fwQXjenLoNR4v0mgTmJMJLedncmkwT0sCXQgbiaCpUCmiAwSkSjgCmBuozJv4K0NICKpeJuKtrsYkzFd0nvrchmf0Y1eSTFMGNCNk4em8q/PtlNeVQvAodIqrnryS1bnFPHApSc0aaq5bEI6h8qq+XhjHgALNufz5qq9/Pj0IQztldDu78e0LdcSgarWALcA84ENwEuquk5E7hOR6U6x+cBBEVkPfAL8UlUPuhWTMceqsqoaHvl4C1fO+oL8ksoG17ILyli3t5jzRqXVn7v1rEwOHK7kha92k1tUweX/WsKGfcU8fs0EvuVntu6pmT1JS4rm5WU5lFfV8ts31jA4NZ6fnGELvh0LXO0jUNV5wLxG5+7xeazAHc6PMV3a397fRHlVLT8/d3jQe+rW1Hp4ZXkOD36wmf0llYQJ/P6tdTxy1fj6MvPXeUcLnTfqyJDLiYO6M3lwdx5bsI2nFu2gqLyaZ743kSlD/C/gFh4mfHtcOk8s3M49b64lu6CcF26aTHSE7f17LLD52cZ0AGVVNfxrwXaeXLSD6Y8sYsO+4hafs2FfMdMeXshdr60hvVssr/xoCrefPYy3V+/jow159eXmr8tlRO9EBvSIb/D8284aRn5JJeXVtbxw0+Rmk0CdyyakU+tRXl6ew3cmpLdY3nQelgiM6QC+3F5AVa2Hm88YQmF5NTMe/ZzZn+/AW2luyuNR7nxlNYVlVTx+zXhe/fFUsgZ250enDWFYWgK/fWMthytryC+pZNmuQ34nYE0e3J2HvjuW138ylePTW97YZWivBCYM6Eb3+Ch+/a3jvvF7Nh2HJQJjOoDPtuQTHRHGT8/M5N3bTuHkoanc+9Z6fvbiKr/JYO7Xe1mzp4i7zz+OaaP71I/AiYoI44FLTyC3uIK/vreRD9bnodqwWaiOiHDxuH5NagqB/PPq8bzxk5PoFh919G/WdDihnkdgjAE+25zPpME9iIkMJyYynKeuz+LvH27hHx9tYcKAblw7ZWB92YrqWv7y3kZG9U3i2+P6NbnX+IxuXD9lIM8s2cnAHvEM6BHHiN6JbRJnmk38OiZZjcCYENtTWM62/FJOzUytPyci3H5WJqcN68n972xo0Gfw1KId7C2q4DcXHNfsmP1fnDecPkkx7DhQynmjetuYfROQJQJjQmzh5nwATh3WcInmsDDhb5ePITk2kp++sJKyqhoOHK7ksU+3cfZxvZg6JNXf7QBIiI7gT5ccT3REGNPH9HU1ftP5WdOQMS7Zuv8wA3vEtbjk88ItB0hLiibTz8Ss1IRo/n75WK59+kvue2s9EeFCeXUtd53fcmft6cN7se7359nmLaZF9i/EGBcs3naAsx9cwPdmL6W4orrZcrUeZdHWA5yS2bPZ5puTM1P58WlDmLM0m+e+3M3VkzKCns1rScAEw/6VGOOCxz7dRmJMBEu2HeTSfy4mu6DMb7nVOYUUlVc3aRZq7GfnDGN8RgoJ0RHcdpYtx2XaliUCY9rYmpwiFm45wC1nDOU/N04kr7iCb//zc1bsPtSk7MItBxCBk4c2394P3k1Znr9pMh/ecVqrdgozJhiWCIxpY48t2EpiTARXTcpg6pBUXvvJScRHR3DFrC94z9kcps5nm/M5vl8y3YMYlx8TGW7DN40rLBEY04a25x/m3bW5XDdlAIkxkYB3Ru7rPzmJUX2TuPn5Fby2wrvbV3FFNSuzCzklM3BtwBi3tThqSER64t1ScqBveVX9vnthGdM5zfpsO1HhYXzvpIbbLnaPj+K/N07ipv8s446Xvqa0qpaeCdHUepRTMwP3DxjjtmCGj74JLAQ+BGrdDceYziu3qIJXV+Rw5cQMUv2048dHR/D0DSdyy/Mr+N0baxnQI474qHDGZXQLQbTGHBFMIohT1V+5HokxndxTi7bjUbjplMHNlomJDOexayZwx0tf89bXezn7uF5ERVgLrQmtYBLB2yLyLWdvAWOMH1vySnj+y91cdEIf+nePC1g2MjyMh747lnH9U5g82JZyNqEXTCK4Dfi1iFQBdTNjVFWT3AvLmM5hTU4R//x0K++tyyUuMpybzxga1PPCw4Tvnzyo5YLGtIMWE4Gqts2yhcZ0cqpKXnElm/NK2JxXwoLN+SzccoDEmAhuPn0o3ztpoI3xN51SUGsNOXsMn+ocfqqqbwf5vGnAw0A48KSqPtDo+g3AX4E9zqlHVPXJYO5tTHvweJQl2w/y4tJsPt20n+KKmvpraUnR3DltONdMHkCSM1TUmM4omOGjDwAnAs85p24TkZNU9e4WnhcOPAqcA+QAS0Vkrqqub1T0RVW9pfWhm45q6c4C1u4p4lBZNYVlVRSVVzN9TF/OOi6t5Sd3EEVl1Tz7xU5eXJZNdkE5STERnD+6D6P6JZHZK5HMtAS/I4OM6YyCqRF8Cxirqh4AEXkGWAkETATARGCrqm53njcHmAE0TgTmGFJT6+HqJ7+kqsYDQHJsJNW1HtbsKeLMEb06zbr4t7+4kk825TN1SA9+ce5wzhvVm5hI26jdHJuCXYY6BShwHre8ualXPyDb5zgHmOSn3KUiciqwGfiZqmY3LiAiM4GZABkZGUG+vAmFnEPlVNV4uH/GKK6aNIDwMGHOV7u567U1rM4pYkz/lCbPeWlZNhv3lXDPRSNDEHFTa/cU8cmmfH553vCgO3+N6cyCGcD8/4CVIjLbqQ0sB/7YRq//FjBQVU8APgCe8VdIVWepapaqZvXsabMwO7IdB0sBGNEniXBn96zzj+9DVERY/dIKviqqa3ng3Y08/fkOlu4saHI9FB77dBuJ0RFcO2VAqEMxpl20mAhU9QVgMvAa8CowRVVfDOLee4D+PsfpHOkUrrv3QVWtdA6fBCYEE7TpuHYe8CaCgT4boifHRnLOcWm8tXof1bWeBuXnrtpLQWkV0RFhPPzhlnaN1Z9t+YeZt3Yf1021DmDTdTSbCERkhPN7PNAHb9NODtDXOdeSpUCmiAwSkSjgCmBuo9fo43M4HdjQuvBNR7PrYBkJ0RGkJjRcTfOS8f0oKK1iwab8+nOqytOf72BE70R+ce5wFm094EqtoKK6lh//dznvrtnXYtnHPt1GdEQY3z/JxvibriNQjeAO5/ff/Pz8b0s3VtUa4BZgPt4P+JdUdZ2I3OcMRwW4VUTWicjXwK3ADUf1Lkyb+XTTfm76z7L6zt7W2nGglIGpcU06hU8d1pMe8VG8vvJIpXDJ9oNszC3h+ycN4prJA0hNiGrzWoGqcvdra3h3bS5Pf74jYNmcQ2W8sXIPV07MsPkApktptrNYVWc6D89X1QrfayIS1KLozrIU8xqdu8fn8d20PPrItJNNuSXc/NwKSqtq2bK/hFF9gx0XcMTOg6WM7tf0eZHhYVw0pi/Pf7WbovJqkmMj+ffnO+keH8X0sX2JiQznR6cN4Q/vbGDpzgJOHNi9Ld4Ssxfv5PWVexjYI47luw5x4HBls8M+n/hsOyKB1woy5lgUTGfx4iDPmU7s4OFKbnxmaf03+c15Ja2+R3Wth5xD5Qzy6R/w9e1x/aiq8TBvzT52HSzlww15XD0po35Y5tWT/NcKKmtqWbuniL2F5Q36GFSV/SUVfLWjgLdX7+XA4coGz/ti+0H+8M4GzhmZxiNXjcej8PGG/X5jyy+pZM7SbC4Zl07flNhWv3djOrNmawQi0hvvENBYERkH1NX1k4DAq2qZTqWqxsOPn1tBfkklz980iStmfcHG3NYnguyCMmo9ysBU/4nghPRkBveM5/UVe9iSd5iIMOGayUdG5sRGNawV9IiP4oWvdvPK8hwOlXmXuRKB1IRokmMj2VtYTlnVkZXRI8OFbx3fh2snD6BvSiw3P7eCAT3iePDyMSRER9AvJZb31+dx+Yn9m8T29Oc7qK718KPTh7T6fRvT2QWaR3Ae3jb7dOBBn/MlwK9djMm0I1XlnjfX8tWOAh6+YiwTBnRnSM8ENh9FIth10LtB+6BU/98TRIRLxvXjf9/fzJo9RVxwfJ8mWy9ePWkAjy/Yxo2zl1JcUUNEmHDOyDSmje5NaWUtucUV7C+uoLCsmlMyUxnYI54BPeJIio1k7qq9vLo8hzdX7SU2MpzwMGHWtVn1O4WdMzKNF77aTVlVDXFRR/7pH66s4b9f7OL80X0Y1EwSM+ZYFqiP4BngGRG5VFVfbceYTDt6bcUe5izN5pYzhjJjbD8AhvdOZOmO1o/e2eEMHR3QTNMQwMVOIiivrvW7+mZsVDh3ThvBUwt38MPT+vKdrHR6JQa3T+/4jG7cOW04b67ay9xVe5l56mCG9kqov37uyDRmL97Jwi0HOG9U7/rzLy7NpqSihpmnWt+A6ZqCmVk8WkRGNT6pqve5EI9pZ2+s2sPgnvHccc6w+nPDeyfy5qq9FFdUt2os/c6DpSRGR9AjwEbs6d3iOH14T6pqPJyQ3nSWMcDlWf25PKtp800w4qIiuHJiBldObDoD/cRB3UmKieCD9Xn1iaCm1sPTi3YwcVB3v7OejekKgkkEh30exwAXYuP9jwlVNR6W7izgihMzCAs7MtxzeJp35fEteSVMGBD86B3v0NH4FtcTmnVt1tEF/A1Fhodx5ohefLQhj5paDxHhYcxbm8uewnLund7ku44xXUYwM4v/5vPzR+B0wOrQx4CVuw9RUe1hypCGu2QN7+1NBK3tMN55sLTZjmJfURFhIdue8dxRvTlUVs3yXYdQVZ74bDuDU+M5a0SvkMRjTEdwNP8b4/B2IJtO7vNtBwkTmmyX2C8lloToiFZ1GFfVeNhzqJxBPTr2gLJTh/UkKjyM99fn8eWOAtbsKeIHpwxuUCMypqsJZj+CNYA6h+FAT8D6B44Bi7ce4Ph+ySTHNuwHEBGGpSW0qkaQfagMjwbuKO4IEqIjmDq0Bx+sz2PngVJ6xEdxyfh+oQ7LmJAKpo/gQp/HNUCes3yE6cRKK2tYlV3ITc2MlBneO5H31uaiqkHtIVC/2FwnGH557sje/Pr1NewuKOP2szNtnwHT5QXTR7AL6IF3U5lLgOPdDsq476sdBdR4lJOGpPq9PjwtkUNl1eSXVPq93ljd0NHOMA7/7OO8/QHREWFcO9mWmjamxUQgIvfg3SegB5AKzBaR37odmGkb76/LZf663CbnP996gKjwMCYM6Ob3ecOcDuNNQS41sfNgKUkxEXSL6/hLN/dKiuHisX354WlDbHE5YwiuaehqYEzdwnPOHsargD+4GZhpGw+8u5Hc4go+/eXpDSZmLd52kPEDUoiN8t8sUjeEdFNuCadktrwZ0K6DZUENHe0oHrpiXKhDMKbDCGbU0F688wfqRNNogxnTMZVV1bDjYCllVbUNFnIrKK1i/b7iZpuFAHokRJOaEM2mIDuMdxwobbAZjTGm8wi0Mc3/icg/gCJgnbNV5b+BtUBhewVoAntx6W625x/2e21jbgmqMLhnPHOWZrPNKbdk20EApg5tPhEADO+dENQqpJU1tewtLO8UHcXGmKYC1QiW4d2f+HW8i8x9AnwK/AZ40/XITIvyiiv41atreGLhdr/XN+wrBuDBy8cSExHGX97bCMDn2w4QHxXOCemB9xsYnpbE5rzDeDxaf27r/sNc9H+LGuwkll3gHTra3GJzxpiOraVF50wHVvfNflV2kd/rG/YVkxgTwZj0ZH502hD+9sFmlu0sYPHWA0wa3IPI8MAtg8N7J1BeXUv2oTIG9IivX6l0zZ4i7nhpFfNvP5W4qAh2HPCuOmpNQ8Z0ToGahl5yfq8RkdWNf4K5uYhME5FNIrJVRO4KUO5SEVERCc0iNJ3U4m0HANiUW0xZVdOpHev3FnNcnyREhBtPGUSvxGjuem0NOw+WMbXRshL+DO+dBBxZauKt1ftYvO0gl01IJ7ugnL+8twmAXQebblhvjOk8Ao0aus35fWGAMs0SkXDgUeAcvJveLxWRuaq6vlG5ROe1vjya1+nKFm87SHJsJEXl1azJKWKSz1IRHo+yMbekfhXPuKgI7jhnGHe9tgaAk1roHwDIdJZw3pxbwtQhPfjD2+s5vl8yf770BBKiI5i9eCfTRvdmx4FSkmMj6RZg1VFjTMfVbI1AVfc5H+azVXVX458g7j0R2Kqq21W1CpiDd1JaY/cDfwYq/FwzzcguKCPnUDnXTx0IwKrshv33uwrKKKuqZWSfpPpzl01IJ7NXAqkJUfXDQwOJj44go3scG/NKeOjDLeQfruT+i0cTHibcOW04Gd3juPOV1WzYV2wdxcZ0YgEbiVW1FvCISOt3Mfduc5ntc5zjnKsnIuOB/qr6zlHcv0ur6x+48IQ+9O8e2yQR1HUUH+eTCCLCw3jq+hP59w0Tg15kbVhaIku2HWT24p1ccWIGY501++OiIvjLZSewu6CMFbsLO/xic8aY5gW7H8EaEfkAKK07qaq3fpMXFpEwvFtg3hBE2ZnATICMjKYbjnRFi7cdIDUhisxeCYzt341lOxvuKLZhXzHhYUJmWkKD8xmt/MAe0TuRDzfk0S0ukjvPG97g2uTBPbh+ygCeWbKrwy82Z4xpXjCJ4DXnx5f6K9jIHsB3m6l0Gk5ESwRGA586s1F7A3NFZLqqLmvwYqqzgFkAWVlZwbz2MU1VWbztIFOGpCIijO2fwltf7yWvuKJ+D+D1e4sZ0jP+Gy+oNrKvt0Zx1/kj/PYB3DltBCWVNUwb3bvJNWNM5xBMIkhR1Yd9T4jIbc0V9rEUyBSRQXgTwBXAVXUXVbUI79pFdff8FPhF4yRgmtp+oJT9JZVMcTqH65prVu4urP9A3rCvmBMHBb+7WHPOHZnGf2+cxElD/Y8yio+O4MHLx37j1zHGhE4wS0xc7+fcDS09yVmq+hZgPt6tLV9S1XUicp+ITG9VlKaBxXUzg50hoKP6JhEZLvX9BIVlVewtqmjQUXy0IsLDODkztdOsIWSMab1mawQiciXeb/CDRGSuz6UkoMD/sxpS1XnAvEbn7mmm7OnB3NPAkm0H6JscwwCnvT8mMpzj+iSxKvsQAOv9dBQbY0xzAjUNLQb24W2++ZvP+RIgqAllpu15PMqSbQc5Y0SvBt/Sx/ZP4dXlOdR6lPV7LREYY4IXaB7BLlX9FDgbWKiqC/AmhnTA2glCZFNeCYfKqpnaaOXQsf1TKK2qZcv+EjbsK6FnYjQ9E22tfWNMy4LpI8GrHPMAABc1SURBVPgMiBGRfsD7wLXAbDeDMs2r6x+Y0miJiLoO46+zC9mwr9hqA8aYoAWTCERVy/BuU/lPVf0OMMrdsExzlmw7wMAecfRLiW1wflBqPMmxkSzdeYgt+0vapKPYGNM1BJUIRGQK3p3K6mYA227fIVBT6+HL7QVNagMAIsKY/im8u2Yf1bXKcX1aXkLCGGMguERwO3A38Loz/HMw3r0JTDvLLa6gpLKGE9JT/F6v6ycArEZgjAlaixPKnE7iBT7H24FvtLyEOTrF5d6lppvbIH6c008QHRHGIFsEzhgTpEDzCB5S1dtF5C38LCmhqjYprJ2VVFQDkBjjPxHU7Tg2vHciES1sOmOMMXUC1QiedX7/b3sEYlpWUuGtESTG+P9r65EQzZj+KfVLTxhjTDACbVW53Pm9oLkypn0VOzWCpGZqBABv/GRqe4VjjDlGBGoaWkOAVUZV9QRXIjLNaqlGANiaQMaYVgvUNFS3ReXNzu+6pqJrCG4ZatPGissD9xEYY8zRCNQ0tAtARM5R1XE+l34lIiuAZjejN+4oqawhJjKMqAjrCDbGtJ1gJ5Sd5HMwNcjnmTZWUlFttQFjTJsLZmOaG4GnffYtLgS+715IpjnF5TUkBegfMMaYoxHMhLLlwJi6RODsLGZCoNhqBMYYFwT99dISQOiVVNQEHDFkjDFHw9r6O5HiimqSYq1GYIxpW64mAhGZJiKbRGSriDQZZSQiPxKRNSKySkQWichIN+Pp7EoqrI/AGNP2Ak0ouyTQE1X1tUDXRSQceBQ4B8gBlorIXFVd71PseVV93Ck/HXgQmBZk7F1Ocbn1ERhj2l6gr5cXBbimQMBEAEwEtjqrlSIic4AZQH0iUNVin/Lx2ES1ZlXVeKis8ViNwBjT5gJNKPveN7x3PyDb5zgHmNS4kIjcDNwBRAFn+ruRiMwEZgJkZGR8w7A6p5ZWHjXGmKMV1NdLEbkA7/aUMXXnVPW+tghAVR8FHhWRq4DfAtf7KTMLmAWQlZXVJWsNxc46Q0mxViMwxrStFjuLReRx4LvATwEBvgMMCOLee4D+PsfpzrnmzAEuDuK+XVJ9jSDaagTGmLYVzKihqap6HXBIVX8PTAGGBfG8pUCmiAwSkSjgCmCubwERyfQ5vADYElzYXU/d7mQ2j8AY09aC+VQpd36XiUhf4CDQp6UnqWqNiNwCzMe72f3Tzp7H9wHLVHUucIuInA1UA4fw0yxkvOpqBDaPwBjT1oJJBG+LSArwV2AF3pE9TwZzc1WdB8xrdO4en8e3BR9q1xbMXgTGGHM0gllr6H7n4asi8jYQY8tNHB2PR6lVJfIo9hMutlFDxhiXBJpQdqaqfuxvYpmItDihzDT15KLtvPBVNp/84vRWP7e4ogYRSIy2GoExpm0F+lQ5DfgY/xPLgplQZhrZkneYHQdKOVxZQ0IrP9BLKqpJiIogLMy2ojTGtK1AE8r+x/n9TSeWGUehs9XkvsJyMtMSW/Xc4vIa6yg2xrgiUNPQHYGeqKoPtn04x7YiJxHsLapodSLw7k5mzULGmLYX6JOl7pNqOHAiR+YAXAR85WZQx6qisiM1gtYqtkRgjHFJoKah3wOIyGfAeFUtcY7vBd5pl+iOMYXlVYC3RtBaJRU19E6KabmgMca0UjDjGNOAKp/jKuecaaWi8qOvEdjuZMYYtwTzyfIf4CsRed05vhh4xr2Qjk0V1bVUVHsAyC1ufY3AdiczxrglmAllfxSR94CTnVPfU9WV7oZ17KmrDQDsbWWNQFWtRmCMcU1QnyyqulxEsnGWoRaRDFXd7Wpkx5hCp6O4d1IM+4oqUFVEgpsTUF5dS61HbVaxMcYVwSxDPV1EtgA7gAXO73fdDuxYU1cjOK5PImVVtfWriQajrmySJQJjjAuC6Sy+H5gMbFbVQcDZwBeuRnUMKizz9reP6JMEwN6i4JuHjuxOZk1Dxpi2F0wiqFbVg0CYiISp6idAlstxHXPqZhWP6O2dnrGvFYmg2BKBMcZFwXyyFIpIAvAZ8JyI7AdK3Q3r2FPsJIKRdTWCwuBHDh3ZptKahowxbS+YGsEMoAz4GfAesA3/C9GZAArLqgkPEwalxhMRJq2qEdTtRZBkNQJjjAuCGT5a9+3fAzwjImHAlcBzbgZ2rCksryIpJoKI8DDSkmLY15oagVObsM5iY4wbmq0RiEiSiNwtIo+IyLnidQuwHbi8/UI8NhSV15ASFwVAn+SYVnYW1+1OZonAGNP2AjUNPYt3wbk1wA+AT4DvABer6oxgbi4i00Rkk4hsFZG7/Fy/Q0TWi8hqEflIRAYcxXvoFArLqkh22vj7pMSyrxXrDRVXVBMRJsREtn5nM2OMaUmgpqHBqno8gIg8CewDMlQ1qE8wEQkHHgXOAXKApSIyV1XX+xRbCWSpapmI/Bj4C/Ddo3gfHV5ReTXdnBpB3+QY5q8LflJZibO8RLAT0IwxpjUCfcWsXxNBVWuBnGCTgGMisFVVt6tqFTAHb8dzPVX9RFXLnMMvgPRW3L9TKSqvJiXOqREkx1BV4+FgaVULz/Ky5SWMMW4KlAjGiEix81MCnFD3WESKg7h3PyDb5zjHOdecG2lmxrKIzBSRZSKyLD8/P4iX7ngKy6pJ8WkaAoLuMC4ut70IjDHuaTYRqGq4qiY5P4mqGuHzOKktgxCRa/BOUvtrM7HMUtUsVc3q2bNnW750u6j1KMUV1STXNw15E0GwHcYlFTU2YsgY4xo3ex/3AP19jtOdcw2IyNnAb4DpqlrpYjwhU1JRjSo+ncXeDWaC3ZfAmoaMMW5yMxEsBTJFZJCIRAFXcGS7SwBEZBzwL7xJYL+LsYRU3YJzdU1DPeKjiIoIC3rkUHFFtdUIjDGucS0RqGoNcAswH9gAvKSq60TkPhGZ7hT7K5AAvCwiq0RkbjO369TqlqCu6ywWEWcuQXCJwFsjsERgjHGHq+0NqjoPmNfo3D0+j8928/U7iroF55J91grqkxwTVNNQrUc5XGlNQ8YY99gMpXZQ3zQUdyQR9E0OblLZYVtwzhjjMksE7aDI2YsgOTaq/lyflBhyiyuo9WjA59oS1MYYt1kiaAd1fQQNm4ZiqfUo+SWBB0rVJQJbedQY4xZLBO2gqLyauKhwoiKO/HH3dYaQtjSX4MgS1NY0ZIxxhyWCdlBYfmRWcZ0+ycHNLq5bgtpGDRlj3GKJoB0Ulh2ZVVynbnZxSxvU1NcIYq1pyBjjDksE7aCovIrkRh/kSbERxEWFt7hl5ZGN661GYIxxhyWCdlBUXk1KbMMaQd2kstziwDWC4vpNaaxGYIxxhyWCdlBYVt1gDkGdPsmxQdUIYiPDiQy3vypjjDvs06UdFJZXk+w3EcQE1UdgtQFjjJssEbisorqWqhpPgzkEdfqkxLK/pJLKmtpmn1/s7E5mjDFusUTgsvoF5xr1EQCMy0hBFT5Yn9fs861GYIxxmyUClxWWe5eX8NdHcGpmTzK6x/Gfxbuafb53dzKrERhj3GOJwGVFfpaXqBMeJlwzOYOvdhawYZ//3T+9u5NZjcAY4x5LBC7ztwS1r8uz+hMdEcZ/lvivFRTbXgTGGJdZInBZUVnTJah9pcRFcfHYfryxck/9ctW+vLuTWY3AGOMeSwQuO7IXQdPO4jrXThlAeXUtryzPaXC+ssY74shGDRlj3ORqIhCRaSKySUS2ishdfq6fKiIrRKRGRC5zM5ZQKSyvIjxMiI8Kb7bM6H7JTBjQjWeX7MTjsz9Bic0qNsa0A9cSgYiEA48C5wMjgStFZGSjYruBG4Dn3Yoj1ArLvCuPikjActdNGcDOg2V8tiW//lzdyqO2BLUxxk1uftWcCGxV1e0AIjIHmAGsryugqjudax4X4wip5mYVN3b+6D7cn7CBpz/fSU2tsmjrgfqkEMzzjTHmaLmZCPoB2T7HOcCko7mRiMwEZgJkZGR888jaUXF5dbMjhnxFRYRx1cT+/OPjrXy2OZ+YyDAmDurBVRMzmDqkRztEaozpqjpF47OqzgJmAWRlZQXe5LeDKSyrJjWh+Y5iXzedOpjEmEhG9UtifEY3YiKb71cwxpi24mYi2AP09zlOd851KYXlVQztlRBU2cSYSG46dbDLERljTENujhpaCmSKyCARiQKuAOa6+HodUlFZcE1DxhgTKq4lAlWtAW4B5gMbgJdUdZ2I3Cci0wFE5EQRyQG+A/xLRNa5FU8o1HqU4ooaSwTGmA7N1T4CVZ0HzGt07h6fx0vxNhkdk4rLA88qNsaYjsBmFruoyBKBMaYTsETgopYWnDPGmI7AEoGLCsu8exEk+9mUxhhjOgpLBC6ypiFjTGdgicBFRdY0ZIzpBCwRuKgwwO5kxhjTUVgicFFhWTUJ0RFEhtsfszGm47JPKBcVBbngnDHGhJIlAhet2VPIgB5xoQ7DGGMCskTgkm35h9mcd5hzR6aFOhRjjAnIEoFL5q/LBeC80b1DHIkxxgRmicAl763NZWz/FPokx4Y6FGOMCcgSgQv2FJazOqeIaVYbMMZ0ApYIXDB/rbdZaNooSwTGmI7PEoEL3luXy4jeiQxMjQ91KMYY0yJLBG0sv6SSpTsLrFnIGNNpWCJoYx+sz0MVSwTGmE7DEkEbe29dLoNS4xmelhjqUIwxJiiuJgIRmSYim0Rkq4jc5ed6tIi86Fz/UkQGuhmP24rKq1m89QDnjeqNiIQ6HGOMCYpriUBEwoFHgfOBkcCVIjKyUbEbgUOqOhT4O/Bnt+JpDx9tyKPGo9YsZIzpVNzcvH4isFVVtwOIyBxgBrDep8wM4F7n8SvAIyIiqqptHcxLS7N5YuH2tr5tA/mHK+mTHMMJ/ZJdfR1jjGlLbiaCfkC2z3EOMKm5MqpaIyJFQA/ggG8hEZkJzATIyMg4qmBS4iLJTEs4qucGKzMtgfNH9yEszJqFjDGdh5uJoM2o6ixgFkBWVtZR1RbOHdWbc22ClzHGNOFmZ/EeoL/Pcbpzzm8ZEYkAkoGDLsZkjDGmETcTwVIgU0QGiUgUcAUwt1GZucD1zuPLgI/d6B8wxhjTPNeahpw2/1uA+UA48LSqrhOR+4BlqjoXeAp4VkS2AgV4k4Uxxph25GofgarOA+Y1OnePz+MK4DtuxmCMMSYwm1lsjDFdnCUCY4zp4iwRGGNMF2eJwBhjujjpbKM1RSQf2NXKp6XSaLZyB2Qxtp3OEKfF2DYsxuANUNWe/i50ukRwNERkmapmhTqOQCzGttMZ4rQY24bF2DasacgYY7o4SwTGGNPFdZVEMCvUAQTBYmw7nSFOi7FtWIxtoEv0ERhjjGleV6kRGGOMaYYlAmOM6eKO+UQgItNEZJOIbBWRu0IdT2Mi0l9EPhGR9SKyTkRuC3VMzRGRcBFZKSJvhzoWf0QkRUReEZGNIrJBRKaEOqbGRORnzt/zWhF5QURiQh0TgIg8LSL7RWStz7nuIvKBiGxxfnfrgDH+1fn7Xi0ir4tISkeL0efaz0VERSQ1FLEFckwnAhEJBx4FzgdGAleKyMjQRtVEDfBzVR0JTAZu7oAx1rkN2BDqIAJ4GHhPVUcAY+hgsYpIP+BWIEtVR+Ndnr2jLL0+G5jW6NxdwEeqmgl85ByH0myaxvgBMFpVTwA2A3e3d1CNzKZpjIhIf+BcYHd7BxSMYzoRABOBraq6XVWrgDnAjBDH1ICq7lPVFc7jErwfXv1CG1VTIpIOXAA8GepY/BGRZOBUvHtcoKpVqloY2qj8igBinR354oC9IY4HAFX9DO+eIL5mAM84j58BLm7XoBrxF6Oqvq+qNc7hF3h3QgyZZv4cAf4O3Al0yNE5x3oi6Adk+xzn0AE/ZOuIyEBgHPBlaCPx6yG8/5A9oQ6kGYOAfODfTvPVkyISH+qgfKnqHuB/8X4r3AcUqer7oY0qoDRV3ec8zgXSQhlMEL4PvBvqIBoTkRnAHlX9OtSxNOdYTwSdhogkAK8Ct6tqcajj8SUiFwL7VXV5qGMJIAIYDzymquOAUkLflNGA08Y+A2/S6gvEi8g1oY0qOM4Wsh3y2yyAiPwGbzPrc6GOxZeIxAG/Bu5pqWwoHeuJYA/Q3+c43TnXoYhIJN4k8JyqvhbqePw4CZguIjvxNq+dKSL/DW1ITeQAOapaV5t6BW9i6EjOBnaoar6qVgOvAVNDHFMgeSLSB8D5vT/E8fglIjcAFwJXd8A9z4fgTfxfO/9/0oEVItI7pFE1cqwngqVApogMEpEovB1zc0McUwMiInjbtTeo6oOhjscfVb1bVdNVdSDeP8OPVbVDfZNV1VwgW0SGO6fOAtaHMCR/dgOTRSTO+Xs/iw7Wod3IXOB65/H1wJshjMUvEZmGt8lyuqqWhTqexlR1jar2UtWBzv+fHGC88++1wzimE4HTiXQLMB/vf7iXVHVdaKNq4iTgWrzfslc5P98KdVCd1E+B50RkNTAW+FOI42nAqa28AqwA1uD9/9chlh8QkReAJcBwEckRkRuBB4BzRGQL3trMAx0wxkeAROAD5//O4x0wxg7Plpgwxpgu7piuERhjjGmZJQJjjOniLBEYY0wXZ4nAGGO6OEsExhjTxVkiMB2SiNQ6wwHXishbrV1VUkTuFZFfOI/vE5Gz2yCmWBFZ4Cxm6CpnJdWfuHj/C0XkPrfubzoXSwSmoypX1bHOKp0FwM1HeyNVvUdVP2yDmL4PvKaqtW1wr5akAH4TgbNg3Tf1DnCRswSC6eIsEZjOYAnOYoEikiAiH4nIChFZ4yzohXPtNyKyWUQWAcN9zs8Wkcucxzvr1oMXkSwR+dR5fJrPhL6VIpLoJ46rcWbXNheHiAx09kJ4wtl34H0RiXWuneism7/KWUd/rXN+lIh85ZxfLSKZeCdvDfEpe7qILBSRuTgzpkXkDqfGtFZEbvd5/Y3Oe94sIs+JyNki8rl49xWYCPVrB32Kd2kG09Wpqv3YT4f7AQ47v8OBl4FpznEEkOQ8TgW2AgJMwDtbNw5Ics7/wik3G7jMebwTSHUeZwGfOo/fAk5yHicAEY3iiQJyfY6bi2Mg3sXPxjrXXgKucR6vBaY4jx8A1jqP/w/vOjl1rxPr3Getz+udjnchvUHOcd37jXfiXYd35dq61z8e7xe95cDTTmwzgDd87nk18H+h/ru2n9D/WI3AdFSxIrKKI8sff+CcF+BPzjISH+KtKaQBpwCvq2qZeldvbe2aUp8DD4rIrUCKHlnjvk4q4Lu/QXNxgHdhuVXO4+XAQKePI1FVlzjnn/e51xLg1yLyK2CAqpY3E+NXqrrDeXwy3vdbqqqH8S5gd4rP669RVQ/eBPGRqirexDHQ53778a6Caro4SwSmoypX1bHAALwfunV9BFcDPYEJzvU8oDXbPdZw5N99/fNU9QHgB3i/jX8uIiMax9PodQLFUelTrhZv7aFZqvo8MN15jXkicmYzRUsD3ceH7+t7fI49jWKJcV7TdHGWCEyHpt4VJW8Ffu50kibj3RuhWkTOwJsoAD4DLnZG9iQCFzVzy514m1UALq07KSJDnG/Rf8a7am2DRKCqh4BwObLHcHNxNPc+CoESEZnknKrfolJEBgPbVfUfePsgTgBK8C6m1pyFzvuNE+8GPN92zrXGMLzNVaaLs0RgOjxVXQmsBq7Eu/FIloisAa4DNjplVgAvAl/j3aVqaTO3+z3wsIgsw/ttvc7tTqfraqAa/ztdvY+3SYbm4mjBjcATTpNXPFDknL8cWOucHw38R1UP4q2ZrBWRvza+kfN+ZwNf4d3R7knnz6k1zsA7esh0cbb6qDFBEpHxwM9U9dqjfH6C056PiNwF9FHV29oyxlbEkgY8r6pnheL1TcfSFuORjekSVHWFiHwiIuF6dHMJLhCRu/H+v9sF3NCmAbZOBvDzEL6+6UCsRmCMMV2c9REYY0wXZ4nAGGO6OEsExhjTxVkiMMaYLs4SgTHGdHH/HwMxglp6cDIpAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "plt.plot(irdf.bins, irdf.rdf)\n", "plt.xlabel('Radius (angstrom)')\n", "plt.ylabel('Radial distribution')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The total number of atom pairs in each distance bin over the trajectory, before it gets normalised over the density, number of frames, and volume of each radial shell, is at `irdf.count`." ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00,\n", " 0.000e+00, 4.000e+00, 2.900e+01, 1.800e+01, 1.600e+01, 2.000e+01,\n", " 8.300e+01, 1.130e+02, 1.080e+02, 1.530e+02, 1.620e+02, 2.090e+02,\n", " 2.430e+02, 2.200e+02, 2.590e+02, 2.690e+02, 3.750e+02, 4.100e+02,\n", " 4.080e+02, 4.760e+02, 4.820e+02, 5.260e+02, 5.630e+02, 6.060e+02,\n", " 6.390e+02, 7.100e+02, 6.780e+02, 7.770e+02, 8.810e+02, 9.440e+02,\n", " 1.003e+03, 1.074e+03, 1.204e+03, 1.237e+03, 1.265e+03, 1.471e+03,\n", " 1.513e+03, 1.526e+03, 1.678e+03, 1.780e+03, 1.782e+03, 2.021e+03,\n", " 1.998e+03, 2.133e+03, 2.309e+03, 2.241e+03, 2.429e+03, 2.535e+03,\n", " 2.713e+03, 2.718e+03, 2.845e+03, 3.021e+03, 3.117e+03, 3.201e+03,\n", " 3.491e+03, 3.673e+03, 3.765e+03, 3.948e+03, 4.096e+03, 4.351e+03,\n", " 4.348e+03, 4.511e+03, 4.748e+03, 4.995e+03, 5.148e+03, 5.503e+03,\n", " 5.600e+03, 5.678e+03, 6.085e+03])" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "irdf.count" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Calculating the average radial distribution function for a group of atoms to itself" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You may want to calculate the average RDF for a group of atoms where atoms overlap; for instance, looking at residue 60 around itself. In this case you should avoid including contributions from atoms interacting with themselves. The `exclusion_block` keyword allows you to mask pairs within the same chunk of atoms. Here you can pass `exclusion_block=(1, 1)` to create chunks of size 1 and avoid computing the RDF to itself." ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "scrolled": false }, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "irdf2 = rdf.InterRDF(res60, res60,\n", " exclusion_block=(1, 1))\n", "irdf2.run()" ] }, { "cell_type": "code", "execution_count": 8, "metadata": { "scrolled": true }, "outputs": [ { "data": { "text/plain": [ "Text(0, 0.5, 'Radial distribution')" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYsAAAEGCAYAAACUzrmNAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3deZzddX3v8ddn5syZNUPWJpCFRBqgiGxGFsUVCYtIaFWqtRKR3rS3qFC9jyu0fRSL2tJrq3Xp1VJIiRZErkJBSsUIuIOQBGRHwhJIIGSZrLNkzvK5f/y+Z+bM5Jw5y5zfnDnJ+/l4nMf5/b7nd37nM1nmc767uTsiIiJjaap3ACIiMvkpWYiISElKFiIiUpKShYiIlKRkISIiJSXqHUAcZs6c6QsXLqx3GCIiDWXt2rXb3H1WodcOyGSxcOFC1qxZU+8wREQaipltKPaamqFERKQkJQsRESlJyUJEREpSshARkZKULEREpCQlCxERKUnJQkRESlKyiNn6LXu4/7nt9Q5DRGRclCxi9i/3PccVtz5a7zBERMZFySJmA6kMewbS9Q5DRGRclCxilso4fYNKFiLS2GJLFmZ2lJk9kvfYbWaXm9l0M1ttZs+G52nhejOzr5rZejN71MxOyrvX8nD9s2a2PK6Y45DOZhlIZclktX2tiDSu2JKFuz/j7ie4+wnAG4E+4DbgCuAed18M3BPOAc4BFofHCuAbAGY2HbgKOAU4Gbgql2AaQSqTBaA/lalzJCIi1ZuoZqgzgOfcfQOwDFgVylcBF4TjZcC3PPIAMNXMDgXOAla7e4+77wBWA2dPUNzjlspENYq+fWqKEpHGNVHJ4oPAd8LxbHd/NRxvBmaH47nAy3nv2RjKipWPYGYrzGyNma3ZunVrLWMfl3SoWfQNqmYhIo0r9mRhZkngfOD/jX7N3R2oSWO+u1/r7kvcfcmsWQX37qiLXM2iV53cItLAJqJmcQ6wzt1fC+evheYlwvOWUL4JmJ/3vnmhrFh5Qxjqs1DNQkQa2EQkiw8x3AQFcAeQG9G0HLg9r/yiMCrqVGBXaK66G1hqZtNCx/bSUNYQ0tlczULJQkQaV6zbqppZJ3Am8Kd5xdcAt5jZJcAG4MJQfhdwLrCeaOTUxQDu3mNmnwMeCtdd7e49ccZdS8M1CzVDiUjjijVZuHsvMGNU2Xai0VGjr3Xg0iL3WQmsjCPGuKVzfRb7VLMQkcalGdwxy9Us+jTPQkQamJJFzIaSheZZiEgDU7KIWa4ZSvMsRKSRKVnELJXNTcpTzUJEGpeSRcxSqlmIyAFAySJG7j602qyShYg0MiWLGOVqFaBmKBFpbEoWMUqH/gpQzUJEGpuSRYxS6fyahZKFiDQuJYsYpfJqFr2aZyEiDUzJIkbpvD4L7ZQnIo1MySJGudnbyeYmrQ0lIg1NySJGuWTR3d6iVWdFpKEpWcQot5fFIe0J+lIZooV1RUQaj5JFjAbTwzULdxhIZUu8Q0RkclKyiNFwzaIF0D7cItK4lCxilA59FrlkoX24RaRRKVnEaHBUslDNQkQaVazJwsymmtn3zOxpM3vKzE4zs+lmttrMng3P08K1ZmZfNbP1ZvaomZ2Ud5/l4fpnzWx5nDHXUm6eRS5ZaBa3iDSquGsWXwF+6O5HA8cDTwFXAPe4+2LgnnAOcA6wODxWAN8AMLPpwFXAKcDJwFW5BDPZ5daGGkoWmmshIg0qtmRhZocAbwOuB3D3QXffCSwDVoXLVgEXhONlwLc88gAw1cwOBc4CVrt7j7vvAFYDZ8cVdy0NhrWhuodqFmqGEpHGFGfNYhGwFfh3M3vYzK4zs05gtru/Gq7ZDMwOx3OBl/PevzGUFSsfwcxWmNkaM1uzdevWGv8o1dmvZqFmKBFpUHEmiwRwEvANdz8R6GW4yQkAj2ap1WSmmrtf6+5L3H3JrFmzanHLcVOfhYgcKOJMFhuBje7+63D+PaLk8VpoXiI8bwmvbwLm571/XigrVj7p5UZDdbepGUpEGltsycLdNwMvm9lRoegM4EngDiA3omk5cHs4vgO4KIyKOhXYFZqr7gaWmtm00LG9NJRNekM1iw7VLESksSVivv8ngBvNLAk8D1xMlKBuMbNLgA3AheHau4BzgfVAX7gWd+8xs88BD4Xrrnb3npjjrolcn0VboolkoknzLESkYcWaLNz9EWBJgZfOKHCtA5cWuc9KYGVto4tfbm2oRHMTHclmDZ0VkYalGdwxyq0NlWxuojOZUDOUiDQsJYsY5daGSjQb7clmdXCLSMNSsojRYOjgTjQZnclm1SxEpGEpWcQoncnS0myYqWYhIo1NySJG6ayTaIr+iNVnISKNTMkiRoPpLIlmAwg1CyULEWlMShYxSmezJJvzaxZqhhKRxqRkEaN0xkfWLDTPQkQalJJFjAYz2eE+i9Zm+lIZormHIiKNRckiRumMk0xEf8QdyQSZrLMvzOoWEWkkShYxSmezJJqiZqiOZDMA/erkFpEGpGQRo8G0k2jO1SyiZKHFBEWkESlZxCgaDZWrWURrNqpmISKNSMkiRtFoqNE1CyULEWk8ShYxikZDjaxZaK6FiDQiJYsYpTPZvNFQUc1Ccy1EpBEpWcQoWhsqqll0toZkkVKyEJHGU3KnPDObBfwPYGH+9e7+sfjCOjBEa0MNz7MA6NunZigRaTzlbKt6O/Bz4MeAvhZXIJ31obWhhpqh1MEtIg2onGTR4e6fqebmZvYisIcoyaTdfYmZTQe+S1RTeRG40N13mJkBXwHOBfqAj7r7unCf5cBfh9t+3t1XVRPPREtnhledVQe3iDSycvos7jSzc8fxGe909xPcfUk4vwK4x90XA/eEc4BzgMXhsQL4BkBILlcBpwAnA1eZ2bRxxDNhUpnh/SySiSYSTaaahYg0pHKSxWVECWPAzPaEx+5xfOYyIFczWAVckFf+LY88AEw1s0OBs4DV7t7j7juA1cDZ4/j8CZMKO+XldGhPCxFpUCWThbtPcfcmd28Lx1PcvbvM+zvwIzNba2YrQtlsd381HG8GZofjucDLee/dGMqKlY9gZivMbI2Zrdm6dWuZ4cUrnXVamof/iDu0p4WINKhy+iwws/OBt4XTn7j7nWXe/3R332RmvwOsNrOn8190dzezmqzZ7e7XAtcCLFmyZFKsA57K2ykPoKO1WTO4RaQhlaxZmNk1RE1RT4bHZWb29+Xc3N03hectwG1EfQ6vheYlwvOWcPkmYH7e2+eFsmLlk14qmx1Vs2jW2lAi0pDK6bM4FzjT3Ve6+0qi/oL3lHqTmXWa2ZTcMbAUeBy4A1geLltONDSXUH6RRU4FdoXmqruBpWY2LXRsLw1lk14646P6LBL0ap6FiDSgspqhgKlATzg+pMz3zAZui0bEkgBucvcfmtlDwC1mdgmwAbgwXH8XUWJaTzR09mIAd+8xs88BD4Xrrnb3HiY5dw8zuEfWLHp6B+sYlYhIdcpJFn8PPGxm9wFG1HdxxdhvAXd/Hji+QPl24IwC5Q5cWuReK4GVZcQ6aaQyUbdJfs2iM5ng5Z6+eoUkIlK1ksnC3b9jZj8B3hSKPuPum2ON6gCQzkbbp+b3WbSrz0JEGlTRPgszOzo8nwQcSjRkdSNwWCiTMeRqFom8ZNGZ1GgoEWlMY9UsPkU0k/qfCrzmwLtiiegAkcrkahbDzVDtyYRqFiLSkIomC3fPTaI7x90H8l8zs7ZYozoApIf6LEbWLAYz2TCzW6vDi0jjKOc31q/KLJM8uZpFbj8LiPosQCvPikjjKVqzMLM5RMtqtJvZiUQjoQC6gY4JiK2hDTdD5dUsWodXnj2kvaUucYmIVGOsPouzgI8SzZj+Ul75HuAvY4zpgJDO7t8MpT0tRKRRjdVnsQpYZWbvc/fvT2BMB4ShZqhRM7hB+3CLSOMpZ1LesWb2+tGF7n51DPEcMApNyhuuWWjJDxFpLOUki715x23AecBT8YRz4EgX6LNQM5SINKpyZnCPmGdhZv9IgyzkV09Dk/KaRu5nAUoWItJ4qhns30HU6S1jKDQpL1ez6FUzlIg0mJI1CzN7jGjGNkAzMAtQf0UJhdaGyiULzeIWkUZTTp/FeXnHaeA1d9dX4xKG14bKW3U2zLNQzUJEGk05fRYbwsKBpxPVMH4BPBx3YI2u0KS81kQTTaahsyLSeMrZVvVvgFXADGAmcIOZ/XXcgTW6QmtDmRkdyYQ6uEWk4ZTTDPVh4PjcYoJhT+5HgM/HGVijK7Q2FET9FppnISKNppzRUK8Qza/IaQU2xRPOgSNVoGYBuWShmoWINJaxNj/6mpl9FdgFPGFmN5jZvwOPAzvL/QAzazazh83sznC+yMx+bWbrzey7ZpYM5a3hfH14fWHePa4M5c+Y2VnV/agTa3g01OiaRUI1CxFpOGM1Q60Jz2uB2/LKf1LhZ1xGNOO7O5z/A/Bld7/ZzL4JXAJ8IzzvcPffNbMPhuv+0MyOAT4IvB44DPixmR3p7pP663mhnfJANQsRaUylFhIcFzObB7wH+ALwKTMzoh32/ihcsgr4LFGyWBaOAb4HfD1cvwy42d33AS+Y2XrgZOD+8cYXp0KT8gA6WhPs6k/VIyQRkaqN1Qx1S3h+zMweHf0o8/7/DPxvIBvOZwA78+ZpbCTaM4Pw/DJAeH1XuH6ovMB78uNdYWZrzGzN1q1bywwvPoXWhgLoaGmmPzRDDaaz/OPdz3DhN+9nX7r82sZ/PryJHb2DtQtWRKSEsZqhLgvP541xTVFmdh6wxd3Xmtk7qrlHJdz9WuBagCVLlniJy2M3vDbU6JpFM737MqzfsofLv/sIj2/aDcD2vYMcNrW95H23793H5d99hL857xg+dvqi2gcuIlLAWM1Qr5pZM3CDu7+zinu/BTjfzM4lGk3VDXwFmGpmiVB7mMfwyKpNwHxgo5klgEOA7XnlOfnvmbRSmSyJJiNqSRvWkWxmy54B3vPVX9DZmuB9J83j++s2lt2PkbtOTVkiMpHGHDobOpGzZnZIpTd29yvdfZ67LyTqoL7X3T8M3Ae8P1y2HLg9HN8Rzgmv3+vuHso/GEZLLQIWAw9WGs9ES2d9vyYogKntSVIZ581HzOCHl7+Vc46dA5S/x0V/KkoWuweULERk4pS7n8VjZrYa6M0Vuvsnq/zMzwA3m9nniZYNuT6UXw98O3Rg9xAlGNz9idB/8iTR2lSXTvaRUBBqFqM6twEufstCTl40nbcunhlmdEfbhVRas9gzoOG3IjJxykkWt4ZHvor6BNz9J4Qht+7+PNFoptHXDAAfKPL+LxCNqGoYqUy2YM1iRlcrbzty1tB5R2tuj4syaxYhWexWM5SITKByksVUd/9KfoGZXVbsYomkM77fsNlCKt09byClmoWITLxylvtYXqDsozWO44CTyviIXfKKaW+pLFnkrlOfhYhMpKI1CzP7ENHkuUVmdkfeS91EfQoyhqgZqnTNIrfHRd++yjq4VbMQkYk0VjPUr4BXiZYlz9+Hew9Q7qS8g1Y6m91vqY9ChpqhUuXVLDQaSkTqYax5FhuADWb2bqDf3bNmdiRwNPDYRAXYqFKZwkNnR8ttiFTuVqu52d97BtK4+37zOERE4lBOn8XPgDYzmwv8CPgIcEOcQR0Iym2Gym2I1Fvm7nn9g9EyIpmsa0FCEZkw5SQLc/c+4A+A/+vuHyBaAVbGkM74fkt9FNOebKY/VVmfBajfQkQmTlnJwsxOI9ox779CWXN8IR0Yis2zKKQz2VxBzWI4QajfQkQmSjm/zS4HrgRuC7OpX0e0ZIeMoZJk0V7BvtwjaxZKFiIyMUpOynP3nwI/zTt/Hqh2qY+DRjrrBZf7KKSzomao7NDx7n41Q4nIxBhrnsU/u/vlZvYDCizv4e7nxxpZgyt3NBREfRbl9j/0D6Zpa2liIJVVM5SITJixahbfDs//OBGBHGjKHQ0FYdny3fvKurY/lWF2dxsbtvexWx3cIjJBxppnsTY8/7TYNVJcOpMta7kPgM5kgt4KFhL8nSmtbNjepz4LEZkwYzVDPcYYq8u6+3GxRHSAqLQZqtxJeX2DGeZP7yDZ3KQ+CxGZMGM1Q+W2U700POeapf6YCpcoPxhV2gxVyaqzHclmprQlitYs/vzGtfzenG4+ccbisuMVERlLqeU+MLMz3f3EvJc+Y2brgCviDq6RVTIaqiOZoD+VIZt1mkpM5OtPZWhvaaa7vaVon8UDz/eUPW9DRKQc5U7Ke0veyZvLfN9BrZJ5FrnFBPvLWEywbzBDW0vxmkU6k2VH3yDbe8vrMBcRKUc5mx9dAqzM24d7J/Cx+EI6MFSULIZ2y8sMLVleTK4ZqrutpeBueTv6UrhDz97ByoMWESmi5G8zd1/r7scDxwPHu/sJ7r6u1PvMrM3MHjSz35jZE2b2t6F8kZn92szWm9l3zSwZylvD+frw+sK8e10Zyp8xs7Oq/WEnUiVrQ3UMbYA0dod1KpMllXHah2oW+1+fq1Fs6x3EXV1LIlIbZTcnufsud99Vwb33Ae8KieYE4GwzOxX4B+DL7v67wA6imgvheUco/3K4DjM7Bvgg0eKFZwP/18wm9dpU7k46W/5oqHK3Vs1tqdqeq1kUaIbaHmoUg+ksvVqVVkRqJLa+B4/sDact4eHAu4DvhfJVwAXheFk4J7x+hkWbNSwDbnb3fe7+ArAeODmuuGshlYm+0Zc9GmqoGWrsmkVueG17cqyaxXDz0/a96rcQkdqItaPazJrN7BFgC7AaeA7Y6e6533IbgbnheC7wMkB4fRcwI7+8wHsmpXQ2Wr+pnJ3yoPyaRa4DPDcaqm8wQyqTHXFNfoLYpn4LEamRsSbl/cFYb3T3W0vd3N0zwAlmNhW4jWiXvViY2QpgBcCCBQvi+piyDNcs4ksWU9qiv7q9A2mmdSaHrtmelyB6epUsRKQ2xhp6894xXnOgZLIYuth9p5ndB5wGTDWzRKg9zAM2hcs2AfOBjWaWAA4BtueV5+S/J/8zrgWuBViyZElde3bT4dt++ZPyymuGyiWTttBnAdGeFiOSRe8+zMBdzVAiUjtjTcq7eDw3NrNZQCokinbgTKJO6/uA9wM3A8uB28Nb7gjn94fX73V3N7M7gJvM7EvAYcBi4MHxxBa3XM2i3LWhyu7gDq93tDSTCjWL0f0W2/cOsmB6Bxu2943ovxARGY9y5llgZu8hGo3Ulitz96tLvO1QYFUYudQE3OLud5rZk8DNZvZ54GHg+nD99cC3zWw90EM0Aoqw4dItwJNAGrg0NG9NWqmKaxYhWZSYdd2fNxoqV3UaPddie+8g86a1s3XPvhFNUiIi41EyWZjZN4EO4J3AdUTf+kt+s3f3R4ETC5Q/T4HRTO4+AHygyL2+AHyh1GdOFulspX0Ww5PyxpJ7vb2lmeYwh2P0kh/b9+7jDfOmMqMrSY9mcYtIjZTz2+zN7n4R0RyIvyXqdzgy3rAaW65mUe7aUM1NRjLRRF+J3fL6R82zgP334d6+d5AZnUmmd7aqGUpEaqacZNEfnvvM7DAgRdTEJEUMN0OVPzK5M9lcshlqIH/obEgW+X0W+9IZ9uxLM7MryczOpJqhRKRmyvltdmcY+vpFYB3wIvCdOINqdOkKJ+VB1BRVdjNUspmu0MGd32eRGyo7o6uV6Z1JLSYoIjVTss/C3T8XDr9vZncCbRUu+3HQGWqGKnM0FOT2tChvBndbopmmJqOrdeQs7lxNYnpnkhldrfSE9aGiifAiItUba1Leu9z93kKT88ysrEl5B6tKJ+VBeRsgDaQytLU0De150d2WGNFnsS3Mq5jZlWRmV5JUxtk9kOaQ9pZKfwQRkRHGqlm8HbiXwpPzKpqUd7DJLfdRSTNUOVur9g1GGx/lTGlrGbGnRa5mMaMzaoaCqGlKyUJExmusSXlXhedxTc47GA2PhqqkgzvB5t0DY17Tn8oMDbMF6G5PjNiHe7jPImqGgmgo7aKZnWXHISJSyFjNUJ8a643u/qXah3NgqHTVWSivZtEfmqFyprS1sGXPcILZ1ruPZHMTXa0JZoSahYbPikgtjNUMNSU8HwW8iWg5DoiapSb1chv1lq6iz6KzjNFQ/YMZ2pPDzVDdbQnWbxnZwT2jK4mZMaMrOVQmIjJeYzVD5Xa2+xlwkrvvCeefBf5rQqJrUMOjoSqrWfSWMRqqo2X4r2z/Pot9Q0liuM9Cw2dFZPzK+eo7G8j/ejoYyqSIaibldYRmqLG2Qu1PZWjLr1m0J9g9kB56T0/vIDM6o76K1kQzU1oT2tNCRGqinIUEvwU8aGa3hfMLGN7RTgqodG0ogM7WBOmsM5jJ0poovGts/2CG2d2tQ+dT2lrIZH2o43vb3kGOmNU19Pr0rqT2tBCRmihnUt4XzOyHwOmh6GJ3fzjesBpbpWtDAUNDYvsHM8WTxejRULn1ofrTtLc0s713uBkKYIZmcYtIjZS1RLm7rzWzlwlLlJvZAnd/KdbIGtjQaKgKZ3BDNJdiakfha6LRUPnzLHJ7WqSY0pZgIJUdGjILML2zlY07+ioNX0RkPyV/m5nZ+Wb2LPAC8NPw/N9xB9bIhnbKS1SwNlRr6d3y+kdNyutuH155dnhC3nDNYmZXUkNnRaQmyvnq+zngVOC37r4IeDfwQKxRNbiq1oZqGXu3PPdc38T+NYvdA+mh5qYRzVBdSXb0DpLN1nWXWRE5AJTz2yzl7tuBJjNrcvf7gCUxx9XQqpmU19EaJYHeIsuUpzJOJuuj5lnk+ixSI5b6yJne2Uo66/vteSEiUqly+ix2mlkX8DPgRjPbAvTGG1ZjS2ezJJqsotVecx3X/UU2QMptfJTfZ9Gdtw93bq+L/JrFzK7hWdxTO4bLRUQqVU7NYhnQB/wF8EPgOQovLihBKuMVjYSCkR3cheSWAulIFu6z2FawZqFZ3CJSGyWThbv3unvW3dPuvgr4OnB2qfeZ2Xwzu8/MnjSzJ8zsslA+3cxWm9mz4XlaKDcz+6qZrTezR83spLx7LQ/XP2tmy6v/cSdGKpOtaCQU5CWLIs1Q/Xm75OW0JppINjexZyBNT+8gHcnmEc1UucSxfa+Gz4rI+BT9jWZm3WZ2pZl93cyWhl/mHweeBy4s495p4NPufgxRB/mlZnYMcAVwj7svBu4J5wDnAIvDYwXwjRDHdOAq4BTgZOCqXIKZrNIZpyVRabIYezTU0MZHecnCzJjSlgh9FiPnWMBwk5RGRInIeI31G+3bRIsIPgb8CXAf8AHgAndfVurG7v6qu68Lx3uAp4C5RM1auRngq4hmhBPKv+WRB4CpZnYocBaw2t173H0HsJoyajb1lMpkK1oXCvJqFqliNYv0iOtyuttb2DOQZnveUh850zrUDCUitTFWB/fr3P0NAGZ2HfAqsMDdx950oQAzWwicCPwamO3ur4aXNjO8ztRc4OW8t20MZcXKR3/GCqIaCQsWLKg0xJpKZbyipT4galJqsjGaoQaj4bjto5LFlLBb3ra9g8yd2jbitWSiie62hBYTFJFxG+s32tB4S3fPABurTBRdwPeBy919d/5rHq2AV5NJAO5+rbsvcfcls2bNqsUtq5bOZisaNgtRk1LHGMuUF+qzgGj4bNRnsW+oQzvfzK5WtqkZSkTGaaxkcbyZ7Q6PPcBxuWMz2z3G+4aYWQtRorgxb8/u10LzEuF5SyjfBMzPe/u8UFasfNJKZbIV7ZKXE+3DXbjPIlfe1rJ/zWJXmGeRv9RHzvTOJD1qhhKRcSr6G83dm929OzymuHsi77i71I0tmmRwPfDUqF317gByI5qWA7fnlV8UOtJPBXaF5qq7gaVmNi10bC8NZZNWKuMV91lALlkUrlnk5lHs12fR1sIrO/tJZ33EUh85M7q0mKCIjF9ZCwlW6S3AR4DHzOyRUPaXwDXALWZ2CbCB4ZFVdwHnAuuJ5nVcDODuPWb2OeChcN3V7t4TY9zjls5kSVY4GgoYuxlqsHAz1JS24feMHg0VlbWydsOOimMREckXW7Jw918Axb5en1HgegcuLXKvlcDK2kUXr/HVLIo0Q+X6LAqMhsoZPRoqKov2tMhmnaYqYhIRgfJmcEuFqu2zaB+rGWowg1k0aipfbjFBKFKz6EySddjZr/WhRKR6ShYxSGedZBXJojOZKD4pLxUtTz56vancYoIQjXwabXqXZnGLyPgpWcQgqlnUtoO7b9ReFjn5NYtpBRYLnNmpWdwiMn5KFjGI+iyqa4bqH2Oexej+Chjus+huSxTsVJ/epVncIjJ+ShYxiEZDVV6z6GxN0FukGWogNXbNolATFAx3emsWt4iMh5JFDKK1oaqoWbQ0M5DKFtzZrm+wSM0i9FkU6twGmNYRvb5NNQsRGQclixhUs58FQGfYLa+/wGKCo/ffzskli0JLfQAkmpuY1tFCj/osRGQclCxikM5mqxoN1R6WKS/UFDVQpM+iKzRDFVrqI2dGV6tmcYvIuChZxKDamkVHqDkU6uQuNhqqucl47/GH8fYjiy+eOLu7lYdf2skrO/srjklEBJQsYlFtn0WuGaq3wDLlxUZDAXztQydy1uvnFL3vp5cexd6BNB/45v28tL2v4rhERJQsYpDOeFVrQ+WaoXIbHeUrNhqqHCctmMZN/+NUegfTfOBff8X6LXuruo+IHLyULGJQzU55AJ253fIqaIYq1xvmHcLNK04lk4U//Nf7+dX6bbzc08e2vfvoG0wTLc0lIlJYnKvOHpTcnXTWq14bCvZvhnJ3+lOZ/ZYnr9TRc7q55U9P5cPX/Zo/uu7XI15bNLOTez71di02KCIFKVnUWDrMkUhWtdxH4Waofeks7tA2zmQB8LpZXfzgE6fzwPPb6RvMMJDKsHbDDm5/5BU27uhnwYyOcX+GiBx4lCxqLJWJ9squpmbRWaRmUWwvi2rN7GrlvOMOGzo/du4h3P7IKzy9ebeShYgUpD6LGktloppFNX0WuWao0UNn+4vsklcrR86eAsBvX9sTy/1FpPEpWdRYOtQsWqragzuq6I3u4M4li9H7b9dKV2uC+dPbeXqzkoWIFKZkUWO5PotqkkVzk9GaaNpvT4taN0MVctTsbp5RshCRImJLFma20sy2mNnjeWXTzWy1mT0bnqeFcgRWvRIAABDqSURBVDOzr5rZejN71MxOynvP8nD9s2a2PK54a2UwneuzqG5UUaE9LYaboeLrYjpqThfPb+tlX7rwEukicnCLs2ZxA3D2qLIrgHvcfTFwTzgHOAdYHB4rgG9AlFyAq4BTgJOBq3IJZrIarllUmyz2X6Z8qGaRjO+v66g53WSyznNbemP7DBFpXLH99nH3nwE9o4qXAavC8Srggrzyb3nkAWCqmR0KnAWsdvced98BrGb/BDSpjKfPAqKaRbEO7rj6LACOnqNObhEpbqL7LGa7+6vheDMwOxzPBV7Ou25jKCtWvh8zW2Fma8xszdatW2sbdQUGc0Nnq1gbCoo0Qw3G3wy1aGYnLc2mTm4RKahuHdwerS9RszUm3P1ad1/i7ktmzSq+Amvc0pnxN0Pt18Gdir+Du6W5iSNmdfHM5t2xfYaINK6JThavheYlwvOWUL4JmJ933bxQVqx80kpnx98MVaxmEWeyADhqzhSNiBKRgiY6WdwB5EY0LQduzyu/KIyKOhXYFZqr7gaWmtm00LG9NJRNWoPpMCmv2ppFa6Jon0WxJcpr5ag5U3hl1wC7B1Kxfo6INJ44h85+B7gfOMrMNprZJcA1wJlm9izw7nAOcBfwPLAe+DfgzwHcvQf4HPBQeFwdyiatcdcsWpoLjoZqbrKqm7bKNdTJrdqFiIwSW4+pu3+oyEtnFLjWgUuL3GclsLKGocVquM+iumTRXmSeRXtLM2bxJouj5nQD8PTmPSxZOD3WzxKRxqIZ3DU2PBqqul/sna1RssjfX6JvMBPrsNmcww5pY0prQv0WIrIfJYsaG2/NoiOZIJP1oaQD0S55cS0imM/MOLJAJ/dgOsuPnthMNqsNkkQOVkoWNTbcZ1H9ch8wcuXZ/nHukleJo+ZM4ZnX9oyo2Xz9vvWs+PZafvTkaxMSg4hMPkoWNZZbG2o8Q2cBevOSRV8qU5ONj8px9Jwp7OpP8drufQC8tL2Pb/70OQBuXbdxQmIQkclHyaLGcmtDVTt0tj23W17eiKiBwQwdE1SzyO1t8XSYnHf1nU/Q0mScf/xh3PfMFnb0Dk5IHCIyuShZ1Nh414aa1dUKwDOb9w6V9acysc+xyMkNn31m8x7uffo1fvzUFj55xmL+7O1HkMo4dz76yoTEISKTi5JFjQ3mOrirXBvq5EXTmTu1nRt/vWGorG8wPWF9FlM7kszubuXRTbv42x88yRGzOrn4LYs45rBujp4zhe+vm9QT6EUkJkoWNZbOjG8/i+Ym449OWcCvntvO+i3RqKSBVHbCahYQzbf4r0dfZcP2Pj57/utJJqJ/Jn9w0lweeXknz23dW+IOInKgUbKosfHslJfzh2+aT0uz8R8PvAQMT8qbKEfN7gLgnGPn8NbFw4syLjthLk0G//mwahciBxslixobHg1V/WzrmV2tnPuGQ/n+2o30DaajZqgJrFmcvngWc6e281fv+b0R5bO72zh98SxuXbdJcy5EDjJKFjWWzmZpbrJxL83xkVMPZ8++NLc9vClqhprAmsXbj5zFL694F/Omdez32h+cOJdNO/t56MVJvUSXiNSYkkWNpTNekwX/3nj4NI6eM4Xrf/ECEP+Ks+Va+vrZdCabuVUd3SIHFSWLGhvMZKseCZXPzPjIaYfz/NZoT+yJrFmMpSOZ4OxjD+Wux15lIJUp/QYROSAoWdRYOuNVj4Qa7YIT5tLVGk3Smyw1C4D3vXEue/alufBf7+fHT742YmkQETkwKVnUWDqbHddIqHydrQned1K05fhkqVkAvPmImXzx/cexo2+QP/nWGs772i/44eOvqtNb5ACmZFFjg2mvWbIAuOjNC5nZleR1szprds9a+MCS+dz76XfwxfcfR99ghj/7j3VcfMND9Gg5EJEDUmybHx2s0tlszZqhAI6Y1cWavz6zZverpZbmJj6wZD6/f+JcbnrwJT5/51Oc+5Wf8/U/OnFo86SBVIabH3yJf/v5C6QyWU593QxOO2IGp71uBofP6Ih9QycRqQ0lixqLRkMdXBW2RHMTF522kJMWTOPSm9bxh9c+wKeXHolhXP+L59m2d5A3LZzGnEPa+dVz27njN9H6UkfO7uJ/vuMI3nvcYSQOsj8zkUajZFFjg5ls1bvkNbpj5x7CnZ84nSu+/xj/54fPAPDWxTP5+Dt/l1NeNwMAd+e5rb386rlt3PjAS/zFd3/Dl1c/y5+9/QiWvn42r+zs58XtfWzY1osDHzp5AbOmtNbxpxIRAGuUkSxmdjbwFaAZuM7dryl27ZIlS3zNmjUTFlu+i//9QbbtHeQHnzi9Lp8/Gbg7dz+xmcOmtnPcvKlFr8tmnR8/9Rr/ct96frNx136vm0FroomPnHo4K952hJKGSMzMbK27Lyn0WkPULMysGfgX4ExgI/CQmd3h7k/WN7L9pbO1mZTXyMyMs489tOR1TU3G0tfP4cxjZvPL9dt5evNu5k/vYOGMThZM72Dz7gG+du+zXP+LF/j2Axv4/RPnMqOzleYmG3o0mdHcBE2h72N3f4qevkF29KbY2T/IgumdnLRgKm88fBqLZnaqj0SkSg1RszCz04DPuvtZ4fxKAHf/+0LXV1uzeHrzbj5x08PjCZWNO/p5w7xDuOVPTxvXfWTYC9t6+fq967nz0VcYzGQZ65+sGUxtb2F6Z5IpbS08v3UvuweijaSmdbQws0u1EzmwveOoWfzVe46p6r0NX7MA5gIv551vBE7Jv8DMVgArABYsWFDVh7QlmlkcVlyt1uLZXZz7htLfqqV8i2Z28k8XHs8/XXg8EDVzZbJOxp1sFrIeHXsWutoSNOf1GWWzznNb97LupR08/NJOdg+k6vVjiEyI2d1tsdy3UWoW7wfOdvc/CecfAU5x948Xur6efRYiIo1qrJpFo4xX3ATMzzufF8pERGQCNEqyeAhYbGaLzCwJfBC4o84xiYgcNBqiz8Ld02b2ceBuoqGzK939iTqHJSJy0GiIZAHg7ncBd9U7DhGRg1GjNEOJiEgdKVmIiEhJShYiIlKSkoWIiJTUEJPyKmVmW4ENFb5tJrAthnBqrRHiVIy1oRhrQzGW73B3n1XohQMyWVTDzNYUm7k4mTRCnIqxNhRjbSjG2lAzlIiIlKRkISIiJSlZDLu23gGUqRHiVIy1oRhrQzHWgPosRESkJNUsRESkJCULEREpSckCMLOzzewZM1tvZlfUO57RzGy+md1nZk+a2RNmdlm9YyrGzJrN7GEzu7PesRRiZlPN7Htm9rSZPRW27J10zOwvwt/142b2HTOLZ/uzymJaaWZbzOzxvLLpZrbazJ4Nz9MmYYxfDH/fj5rZbWY2dbLFmPfap83MzWxmPWIby0GfLMysGfgX4BzgGOBDZlbdBrbxSQOfdvdjgFOBSydhjDmXAU/VO4gxfAX4obsfDRzPJIzVzOYCnwSWuPuxRMvyf7C+UQFwA3D2qLIrgHvcfTFwTzivpxvYP8bVwLHufhzwW+DKiQ5qlBvYP0bMbD6wFHhpogMqx0GfLICTgfXu/ry7DwI3A8vqHNMI7v6qu68Lx3uIfsHNrW9U+zOzecB7gOvqHUshZnYI8DbgegB3H3T3nfWNqqgE0G5mCaADeKXO8eDuPwN6RhUvA1aF41XABRMa1CiFYnT3H7l7Opw+QLTTZt0U+XME+DLwv4FJOepIySL6pfty3vlGJuEv4hwzWwicCPy6vpEU9M9E/9iz9Q6kiEXAVuDfQ1PZdWbWWe+gRnP3TcA/En3DfBXY5e4/qm9URc1291fD8WZgdj2DKcPHgP+udxCjmdkyYJO7/6besRSjZNFAzKwL+D5wubvvrnc8+czsPGCLu6+tdyxjSAAnAd9w9xOBXurfbLKf0O6/jCi5HQZ0mtkf1zeq0jwahz8pvxUDmNlfETXp3ljvWPKZWQfwl8Df1DuWsShZwCZgft75vFA2qZhZC1GiuNHdb613PAW8BTjfzF4kasp7l5n9R31D2s9GYKO752pl3yNKHpPNu4EX3H2ru6eAW4E31zmmYl4zs0MBwvOWOsdTkJl9FDgP+LBPvsllRxB9MfhN+P8zD1hnZnPqGtUoShbwELDYzBaZWZKoI/GOOsc0gpkZUTv7U+7+pXrHU4i7X+nu89x9IdGf4b3uPqm+Dbv7ZuBlMzsqFJ0BPFnHkIp5CTjVzDrC3/0ZTMKO+OAOYHk4Xg7cXsdYCjKzs4maR8939756xzOauz/m7r/j7gvD/5+NwEnh3+ukcdAni9Dx9XHgbqL/kLe4+xP1jWo/bwE+QvRt/ZHwOLfeQTWoTwA3mtmjwAnA39U5nv2Ems/3gHXAY0T/T+u+HISZfQe4HzjKzDaa2SXANcCZZvYsUY3omkkY49eBKcDq8H/nm5MwxklPy32IiEhJB33NQkRESlOyEBGRkpQsRESkJCULEREpSclCRERKUrKQhmVmmTAU8nEz+0Glq4ma2WfN7H+F46vN7N01iKndzH4aFqiMVVhB989jvP95ZnZ1XPeXxqJkIY2s391PCCuz9gCXVnsjd/8bd/9xDWL6GHCru2dqcK9SpgIFk0VYgHC8/gt4b1iOQg5yShZyoLifsACkmXWZ2T1mts7MHguLtBFe+ysz+62Z/QI4Kq/8BjN7fzh+MbefgJktMbOfhOO3502KfNjMphSI48OEWczF4jCzhWEvjX8Le1b8yMzaw2tvCvsuPBL2YXg8lL/ezB4M5Y+a2WKiCXBH5F37DjP7uZndQZiZbmafCjWvx83s8rzPfzr8zL81sxvN7N1m9kuL9qU4GYbWevoJ0TIZcrBzdz30aMgHsDc8NwP/Dzg7nCeA7nA8E1gPGPBGohnRHUB3KP9f4bobgPeH4xeBmeF4CfCTcPwD4C3huAtIjIonCWzOOy8Wx0KiBe1OCK/dAvxxOH4cOC0cXwM8Ho6/RrSuUe5z2sN9Hs/7vHcQLY64KJznft7OEO8TRCsW5z7/DURfGNcCK0Nsy4D/zLvnh4Gv1fvvWo/6P1SzkEbWbmaPMLw09upQbsDfhSU9fkxU45gNvBW4zd37PFq1t9I1wH4JfMnMPglM9eE9EnJmAvn7YxSLA6KFAh8Jx2uBhaHPZYq73x/Kb8q71/3AX5rZZ4DD3b2/SIwPuvsL4fh0op+31933Ei1I+Na8z3/M3bNESeQed3ei5LIw735biFa+lYOckoU0sn53PwE4nOgXc67P4sPALOCN4fXXgEq2JU0z/H9j6H3ufg3wJ0Tf6n9pZkePjmfU54wVx7686zJEtZCi3P0m4PzwGXeZ2buKXNo71n3y5H9+Nu88OyqWtvCZcpBTspCG59FKop8EPh06dg8h2lsjZWbvJEomAD8DLggjlqYA7y1yyxeJmnAA3pcrNLMjwrfxfyBarXhEsnD3HUCzDe+XXSyOYj/HTmCPmZ0Sioa2UjWz1wHPu/tXifpEjgP2EC2QV8zPw8/bYdEmT78fyipxJFHTmBzklCzkgODuDwOPAh8i2txmiZk9BlwEPB2uWQd8F/gN0W5pDxW53d8CXzGzNUTf+nMuDx3FjwIpCu+49iOi5h+KxVHCJcC/hea1TmBXKL8QeDyUHwt8y923E9VwHjezL46+Ufh5bwAeJNpZ8brw51SJdxKNipKDnFadFakhMzsJ+At3/0iV7+8K/QuY2RXAoe5+WS1jrCCW2cBN7n5GPT5fJpdajMUWkcDd15nZfWbW7NXNtXiPmV1J9H9zA/DRmgZYmQXAp+v4+TKJqGYhIiIlqc9CRERKUrIQEZGSlCxERKQkJQsRESlJyUJEREr6/7jGq4Q7Mm0xAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "plt.plot(irdf2.bins, irdf2.rdf)\n", "plt.xlabel('Radius (angstrom)')\n", "plt.ylabel('Radial distribution')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Similarly, you can apply this to residues. " ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "There are 11 THR residues\n", "THR has 14 atoms\n" ] } ], "source": [ "thr = u.select_atoms('resname THR')\n", "print('There are {} THR residues'.format(len(thr.residues)))\n", "print('THR has {} atoms'.format(len(thr.residues[0].atoms)))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The code below calculates the RDF only using contributions from pairs of atoms where the two atoms are *not* in the same threonine residue." ] }, { "cell_type": "code", "execution_count": 10, "metadata": { "scrolled": true }, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ "irdf3 = rdf.InterRDF(thr, thr,\n", " exclusion_block=(14, 14))\n", "irdf3.run()" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "Text(0, 0.5, 'Radial distribution')" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX4AAAEGCAYAAABiq/5QAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3deXxcdbn48c8zk0z2pW3SNum+l0IXoOybZRUooAICIqLIBRXBBa+K/MR7vder94I7KgIiKIvKJiDIIkvLptCF0r2FrmnSJs2+z/b9/THnpJN0ZjJJzsycNM/79ZpXJieTOd+0yTPfec7zfb5ijEEppdTI4cn0AJRSSqWXBn6llBphNPArpdQIo4FfKaVGGA38Sik1wmRlegDJKCsrM1OnTs30MJRSalhZuXLlfmNMed/jwyLwT506lRUrVmR6GEopNayIyM5YxzXVo5RSI4wGfqWUGmE08Cul1AijgV8ppUYYDfxKKTXCaOBXSqkRRgO/UkqNMBr4Xa7TH+LxlVVo+2yllFM08Lvcixv2cvOja9i+vz3TQ1FKHSI08Ltca1cQgA5/KMMjUUodKlIW+EXkPhGpFZF1UcdGi8hLIrLV+jgqVec/VHRaAb8zoIFfKeWMVM747wc+2ufYt4GXjTGzgJetz1UCdsDv1Bm/UsohKQv8xpjlQEOfwxcBD1j3HwA+lqrzHyo6dMavlHJYunP844wxNdb9vcC4NJ9/2OmyAn6XBn6llEMydnHXROoT49Yoish1IrJCRFbU1dWlcWTu0uGPXNzVVI9SyinpDvz7RKQCwPpYG++Bxpi7jTGLjTGLy8sP2kdgxOgMhK2PGviVUs5Id+B/Grjaun818FSazz/sdNozfg38SimHpLKc8xHgbWCOiFSJyOeBHwFnichW4Ezrc5WAHfC7NNWjlHJIyrZeNMZcEedLZ6TqnIcirepRSjlNV+66nC7gUko5TQO/yx1YwBXO8EiUUocKDfwuZ8/4tY5fKeUUDfwup6kepZTTNPC7nPbqUUo5TQO/iwVCYYLhyOJmnfErpZyigd/Fonvwa45fKeUUDfwuFh3sNfArpZyigd/F7Bl/TpZHUz1KKcdo4Hcx+4Lu6AKfXtxVSjlGA7+LdQYiDdpGF/joCugCLqWUMzTwu5i9Wnd0gQ9/KEwwpMFfKTV0GvhdzN6EZXSBD4CuoAZ+pdTQaeB3MfuC7qj8SODXPL9Sygka+F3MDvRj7Bm/VvYopRyggd/Femb8VuDXkk6llBM08LtYR1Q5J2iqRynlDA38LtYVCCECJXnZgM74lVLO0MDvYh3+EPnZXvJ8XkADv1LKGRr4XawzECLP5yUvOxL4dcN1pZQTNPC7WKe/d+DXGb9Sygka+F2s0x8iT1M9SimHaeB3sY5AiDxfFrn2jF9TPUopB2jgd7Euf4i8bM+BHL/O+JVSDtDA72IdgSD5viyyvYLXI5rqUUo5QgO/i9k5fhEhL9urrZmVUo7QwO9idlUPQG62V2f8SilHaOB3sc5AqCe/n+fzaB2/UsoRGvhdrMMfIt+a8efpjF8p5RAN/C4VDhu6g+GeUk4N/Eopp2Qk8IvI10RkvYisE5FHRCQ3E+NwMzvI50fn+DXVo5RyQNoDv4hMAG4CFhtjjgC8wOXpHofb2YHfvrib5/NqHb9SyhGZSvVkAXkikgXkA9UZGodr2bP7PE31KKUclvbAb4zZA9wB7AJqgGZjzIt9Hyci14nIChFZUVdXl+5hZtxBM34N/Eoph2Qi1TMKuAiYBlQCBSLy6b6PM8bcbYxZbIxZXF5enu5hJiUcNlx93zu8sXW/489t775l5/hzsr10+nUBl1Jq6DKR6jkT2G6MqTPGBIAngBMzMI4h6wiEWLaljpU7Gx1/bjvVE13Vk64cf1VjB8u3jLx3WUqNFJkI/LuA40UkX0QEOAPYmIFxDJk/GJmB+0POB+SunqqeLCCygKszEMIY4/i5+vr1ax9y/R9XpuVcSqn0y0SO/1/AY8AqYK01hrvTPQ4n9AT+oPMpmI4YF3dDYUMglPpgvLO+nc5AiNbuYMrPpZRKv6xMnNQY8z3ge5k4t5MCobD10flgHKuO3z7uy0rt6/XO+g4Aalu6Kc7NTum5lFLppyt3h6Dbmul3p2DG3+mPzLZ7cvy+9PTk9wfDVDd1AlDX2p3ScymlMkMD/xCkMtXTd8afl6ZduPY0dRK23sDUtWngV+pQpIF/CPyhcK+PTuqIUdUDqd93d2d9e899nfErdWjSwD8EB2b8zgfjzkCInCwPXo8AkJumDdd3N3T03K9t7UrpuZRSmaGBfwhSmuqJ2oQFSNu+uzvrO8jN9lBZkjugGf+u+o7+H6SUcgUN/ENg1++nItVjb7toS1vgb+hg8uh8youTD/yb9rZw6u2vskwXfSk1LPRbziki5cC/AVOjH2+MuSZ1wxoeUlrHH+gz47dTPSlu27CrPhL4QahqTG4Wv2N/5HGvbqrltNnubK+hlDogmTr+p4DXgX8A2iUsSncKA39XnBl/KnP8xhh2NXRw0swyuoIhVu9KrhXFfqv65/WtOuNXajhIJvDnG2O+lfKRDEP+FNbxR2+7CL0XcKVKXVs3nYEQU8bk09Dup6HDTyAUJtubOCNY3+YH4MO6dmqaO6koyUvZGJVSQ5dMjv9vInJeykcyDKWynLMzEOoJ9hC1gCuFdfz2BdrJY/IpL8rBGGho9/f7ffuj6v1T0alUKeWsZAL/V4gE/y4RabVuLake2HCQ6qqeXjN+q01DKmf8dquGKaMjgR+Sq+Xf39bN9PICygp9vPmBBn6l3K7fVI8xpigdAxmO7IAfSNGMPzrHn+X14PN6Uhr4dzV0IAITRuXR3BkA7Fr+koTft7+tm/LCHOZPKOGND+oxxhBpvKqUcqOkyjlF5EIRucO6LU31oIYLO+Cnqjtnnq/363JutielLRt2NXRQWZJHTpZ3QDP++jY/ZUU5nDyzjP1t3Wza25qyMSqlhq7fwC8iPyKS7tlg3b4iIj9M9cCGg1Smerr6zPghcoE3lXX8O+vbrVJOKCtMPvDXtXVTVuDj5FllAJruUcrlkpnxnwecZYy5zxhzH/BR4PzUDmt46E7RxV1jDB3+YK8cP0Qu8KY61TNlTCTw52Z7Kc7N6jfwdwVCtHYFKSvMoaIkjxnlBbyuF3iVcrVkV+6WRt1PnPAdQQ7k+A3hsHM9+f2hMGFDrwVcYG24nqJUT3t3kP1tfiZZM36AscW51PYT+O2qnzIrNXTyzDLe2d5Adwr6FymlnJFM4P8hsFpE7heRB4CVwA9SO6zhITrF4+Ssv7PP7lu23OzUzfh3Wc3Z7Bk/QHlhTr8zfruU004NnTyrnM5AiFU7m1IyTqXU0PUb+I0xjwDHE9kU/XHgBGPMn1M9sOEgZYHfCu6xZvypyvEfKOUs6DlWXpTTb09+O/CPKfQBcNz00Xg9onl+pVwsbuAXkbnWx6OACqDKulVax0a86GDv5AVeuxd/OnP8uxoiffgnR8/4i5KY8bdGUj3l1oy/ODebRZNKeV0Dv1KulaiO/+vAdcCPY3zNAKenZETDSK8Zv4OBv7PPJiy2VOb4dzV0UJKXTUnegT12xxbl0OEP0dYdpDAn9q/K/vbeM36Ak2aWcecrW2nuCFCSr3v2KuU2cWf8xpjrrLvnGmOWRN+IVPqMeCkL/IHYM/5IOWdqunPurO/old8Hkqrl39/qJ9/nJT9qzcEps8oIG3h7m876lXKjZC7uvpXksRGnV6onDRd383ypW7m7q6Gjp4bfllTgb+vuubBrWzCxhCyP8H5Vs/MDVUoNWdxUj4iMByYAeSJyJGCvwS8G8uN930jSnaIZv53jT1c5ZzAUZk9jJ0sXVPQ6nkzgr2/vpiwqzQOQk+Vl5thCNtRoSyel3ChRjv8c4LPAROAnUcdbge+kcEzDRqqqeuzKnYNm/NleuoIhx3vhVDd1EQybXhU9AGOLcoHEe+/ub/UflCICmFdZrJ06lXKpuIHfGPMA8ICIXGyMeTyNYxo2/MEw+T4vHf5Qiqp6+vTq8XkxJvJOo++F36Gwa/gn9Un1lOZlk+WRflM9R00ZddDxeRXFPLFqD3Wt3T3vHJRS7pDMRixHiMjhfQ8aY76fgvEMK4FQmMKcLMcDf2eCGT9E3hE4Gfh3WqWcfWfuHo9QlmARVzAUpqHDT3mfVA9EZvwAG2taKC/S7RiVcpNkLu62Ae3WLQScS2T/3RHPHwpTmBt57XS2nDMIxM7xg/M9+XfVd+DL8jC+OPegryVaxNXYEcCYA+0aos2riAR+zfMr5T7J9OPvVccvIncAL6RsRMOIPxhmrBX0nF656/UI2d7eefwDG647G/irGjuZWJqHx3PwdYOxRTnUNMfO8fdt1xCtNN/HhNI8NlRr4FfKbZJt0hYtn8gF3xHPH0zNjL/DHyI/23vQBdxU7bvb2OFndMHB6RpIPOPvadcQ53sPqyhmo874lXKdZPrxrxWR963bemAz8LOhnFRESkXkMRHZJCIbReSEoTxfpviD4Z4VrU4G/q5AiFzfwTn86By/k1q6AhTnxV5hW16UQ31bN6EY3UftTdZjpXogkuf/sK4tpXsIKKUGLpmLu9E7bgWBfcaY4BDP+3PgeWPMJSLiY5iuC+gOhSnMye6575S+++3aDqR6nF2929wZYNbY2DtslhflELY2Xe9bnZMo1QORPH/YwOa9rSycVBrzMUqp9EumO+dOYAxwEfAJYP5QTigiJcCpwO+s5/cbY4ZdD19jDP5gmKIUpXr6VvQA5GalJtXT0hmkODf2HMC+hhGrlr+urRuf1xP3ew+v1Au8SrlRMqme24AHiAT/MuB+Efl/QzjnNKAO+L2IrBaRe0WkoO+DROQ6EVkhIivq6uqGcLrUCIQiqY9UpHo6A6GDKnog0rLB/rpTwmFDS1egV3O2aIlW7+5v9TOm0Bd3MdnEUXkU5WTpBV6lXCaZi7tXAscYY75njPkekd78Vw3hnFnAUcBvjDFHEikT/XbfBxlj7jbGLDbGLC4vd18duF3FU5CKwB9vxm/n+B2s6mnzBzGG+Dn+wkiJZ6zAH2nXEH9xlohwWGWxzviVcplkAn81EF3gnQPsGcI5q4AqY8y/rM8fI/JCMKzYgT4v24PXI/hDzgXjzkCcHH8KqnqaOwJAgsBvz/hjVPZEGrTFruixzbMqe5zcmlIpNTSJmrT9kkjf/WZgvYi8ZH1+FvDOYE9ojNkrIrtFZI4xZjNwBrBhsM+XKXbg92V58Xk9PakfJ3T6Y6/M7bm462Tg77QCf27swJ/n81KUk0VtS+xUz9zxxQmff15lMR3+EDsbOphWdlBGTymVAYmqelZYH1cCT0Ydf82B894IPGRV9GwDPufAc6bVgcDvwZflcTzHH2vG33Nx18FUT0tXJPDHy/FD7Fp+Y0y/qR6IWsFb3aKBXymX6K9JW0oYY94DFqfq+dPBTu1ke4Vsr6dXi+ahilfV4/EIOVkeR+viWzr7D/xlMbZgbOkMEgiZflM9s8YVkuURNtQ0c36fts9KqcxIlOr5izHmkyKylkiKpxdjzIKUjszl7ECfk+UhJwUz/jxf7P8ap/fdbemMLMkozov/5q+8KIeNfSpz7HcA/XXe7OnNr5U9SrlGolTPV6yPSxM8ZsSyc/o9qR6HFnCFwpH1AbFm/OD8ZizNScz4xxblsKzPjP9Au4b+Wy7PqyzmTd18XSnXSLTnbo2IeIH7jTE7+97SOEZX6snxeyMXd/1BZ4JxvP12bXnZDs/4uwJ4BArivMOAyIYsbd1Bmjr8PccOtGtInOqBSJ5/X0t3z4uFUiqzEpZzGmNCQNhabauipOribofVkjlWrx6wN1x3dsZfnJcdszOn7aSZYwB4bu3enmP9tWuIFt2bXymVecn2418rIr8TkV/Yt1QPzO3si7tOp3q6rD48+fFSPT4vXQHnric0dwbilnLa5k8oYdbYQh5fVdVzbH9bNx6BUfnJzfgB1u3RwK+UGyQT+J8AvgssJ1LauZIDpZ4j1oFUj8dK9Tg04w/E3oTF5niqpzN+uwabiPCJoyaycmcjO/ZHduva39bN6AIf3gTvFGyl+T5mjS3kjQ/c13pDqZEomcBfaox5IPoGHLzJ6gjTnaJUj33hNt7F3dwUXNxNVNFj+/iRE/AIPGHN+ve3+ZNK89hOnzuWd7Y30GqtG1BqOGntCvD9Zzaw29qferhLJvBfHePYZx0ex7Djjyrn9GU5V8ffE/jjzfh9zub4W7qC/c74AcaX5HLSzDKeWL2HcNhY7RqSD/xL5o4lEDK8sVWre9TwYozh1ifXcd+b27njxc2ZHo4j4gZ+EblCRJ4BponI01G314CGtI3QpeycvtM5/ngbrdvysj2Ot2zoL8dvu/ioiVQ1dvLOjgb2t3Uzpp/FW9GOnjKK4twsXtlUO9ihKpURj66s4uk11Uwancff3q+hqnH4z/oTzfjfAn4MbLI+2revA+ekfmjuFp3jz/F6CDgU+Dv8aS7nTCLHbzvn8PEU+Lw8saqK+gGmerK9Hk6dXc6rm+u0YZsaNj6obeN7T63nhOljePja4xHg92/uyPSwhixRHf9OY8xrwJnA68aYZUANkf12+7+id4hLVTmnHdRjNWmzjzuV4+8KhOgOhuN25uwrz+flvPkVPLOmhg5/aECBHyJ5/v1t3ayrbh7McJVKq65AiC8/vIo8n5efXb6ISaPzWbqggj+9s6tn4eNwlUyOfzmQKyITgBeJ9OK/P5WDGg5SFvj7mfHnZnvpDoYdmTXbfXqSDfwAFx89sefFqb8+PX2dNrscETTdo4aFHz63kU17W7nj0gWMK450pr/u1Bm0+0M89K/hvYY1mcAvxpgOItsu/toYcylweGqH5X52Tj/LE2nS5lTgb7cWcOUn6NUD0OXASuFkOnP2dezU0UwozQOSW7wVbUxhDosmlfKqBn7lUs0dAR7+1y4u+c1bPPD2Tj5/8jROnzuu5+vzKos5ZVYZv39zB90OrdbPhKQCv4icQGQnrmetY7GnoyOIPxjGl+VBRBy9uFvT1EVJXnbCOn5wpjXzgV78/Zdz2jwe4eKjJgADD/wAp88Zy5qq5pg7eimVKXWt3dzw8CqO+Z9/8J0n19LUGeCWc+fyrY/OPeix1506nbrWbp56rzoDI3VGMn/xXwVuAZ40xqwXkenAq6kdlvv5Q2FyvJHXTXsjlnDYJGx9kIw9TZ09M+pYnNyFy+7MOZAZP8A1J08jz5fV04phIJbMHcuPX9rCq5tr+eTiSQP+fqVS4ft/28BLG/bxqWMnc/FREzliQnHcvaRPnlnGYRXF3LN8G5ccNXHIf/OZ0O+M3xizzBhzoTHmf63Ptxljbkr90NzNnvEDPR+dmPVXNXYwcVT8wG/38HGilr95EDl+iKzE/eJHZiS1arevwyuLGVeco+ke5RqrdzXyzJpqrj91Ov9x4eHMn1gSN+hDZCX7dadOY2ttG//YuC+NI3VOojr+n1kfn+lTx/+0iDydviG6U3Tgz3Eo8BtjqGrsZOKo/LiPOZDqGfqLzGBy/EMlIiyZM5bXt+53dA8DpQbDGMMPnt1IWWEO1582I+nvW7qgkhnlBXz3qXXUx+g62x0Mcc397/KNR9dgjPvKlxPN+P9ofbyD3nX89m1E84dizPiHGMgaOwJ0+EMJZ/x24Le7eA5Fz0brSS7gcsqSuWNp6w6yYseIXweo0uRf2+r51asfHLTe5vl1e1mxs5Gbz55NYU7y17qyvR5+ccWRNHYE+PfH3u8V3I0x3PLEWl7ZVMtjK6v4zbIPHfs5nJKojn+l9XFZrFv6huhO/mAYX1SO3z42FPaKwESB3/7aegd2tGrpCpCX7e154UqXk2eW4fN6eG2LNm1TqVfX2s0XHlzJ7S9s5rO/f6dnwuMPhvnR85uYPa6QS4+eOODnPbyyhFvPO4xXNtVyX9Sirt8s+5AnVu3hq2fO4oKFldzxwmZe3+qu3/VEqZ61IvJ+vFs6B+lGMXP8Qw78nQBMSBD4p5YVMHd8Ec+urRnSuSD5Bm1OK8jJYnp5Advq2tJ+bjWyGGP4f39dS7s/xNfPms072xv4+K/fZPv+dv7w9g521nfwnfMOI8s7uMnPZ06YwlnzxvGjv29kbVUzz6/by/89v5kLFlbylTNm8b8Xz2fW2CJuemS1q1o9JPpplwIXAM9btyut29+B51I/NHeLmeoZYo5/jxX4E+X4AS5YWMnKnY1UN3UO6XzNA2jX4LSKklxqmruSfvzGmpZhv2hGpd/Ta6p5Yf0+bj5rNjedMYuH/+14mjoDfOxXb/KLl7dyyqwyPjJn7KCfX0S4/ZIFlBXm8IUHV/K1P7/Hokml3H7JAkSEfF8Wd111NMGQ4YsPrnK0weJQ9NeyYSdwljHmm8aYtdbtW8DZ6RuiO3WnKNVTlJvVbzA+b34FAM8Ncdbf0hlMe37fVlGaN6DA/4uXt3Lrk+vYM8QXO3Vo2N3QwTvbG3hvdxPrq5vZuq/1oLUt+1q6uO2p9Rw1uZRrT5kOwDFTR/PUDScxrjiHtu4g3znvsCGPpTTfx88vP5Ka5k5G5Wdz92eO7tVyZVpZAT+9bBFr9zRzyxNrXdGrKpn3+SIiJxlj3rQ+OZHkFn4d0vzBMEXWwienZvz9VfTYppUVcHhlMc+uren5hR6M5s4AFSW5g/7+oagozqWh3U9XIBS3L5EtGArzhrVZ+/Pr9vL5k6elY4jKpQKhMBfe+QaNHb375RTmZHHBwkquOHYS8yeU8J0n1tIdDHHHpQt7lR5PGp3PX284ieqmTmaOLXJkTMdOG82frjuBytJcxhYd/Dd15rxx3HzWbH780hZE4PZLFg6qHNopyQT+zwP3Re272wRck7ohDQ/+YLinjNPJHP/kMf0HfoDzF1Twf89vtur+k/uevlq6Aswd78wv/kBVWIvU9jZ3MbWsIOFj39vdRGtXkCyP8Pe1NRr4R7gVOxpp7AjwzY/O4bCKYvzBMF2BEMu37OfJ1VU88s4upozJZ2d9B7ctncf08sKDniPfl+VY0LcdO210wq/feMYsDPCTl7YQCBl+8smFZA/y2sJQ9Rv4reqehXbgN8Zoa0V65/hzHAj8kRr+Dk6YMSapx58/PxL4n1tbw3WnJl9/HM3eaD0T7HcaNUkE/mVb6vAIfPbEqdz7xnb2tXT1NM1SI8+rm2vJ9gpXnzCVgqgSzIsWTeC2C+bx9Jpq/vLubs48bCyfPXFq5gYaw01nzMKX5eFHf99EIBjmF1ccmfaqOhhAysYY06xB/wB/MNzzap3tQI6/uTNAez81/NGmjClg/oQSnn1/cHn+cNjQ1h10QeDvP2e/fEsdR04exeXHRlo8vLB+b0rHptzt1U21HDdtTK+gbyvJy+aq46fwzI0nc+/Vx7iyncIXTpvBd5fO4/n1e/nCgytp7x76mpyBGvG5+sHqVcfvQI6/KsmKnmhLF1Swpqp5UPuAtnYFMWZgDdqcVFESeYHr7wJvQ7uf9/c0c+qscmaOLWLW2MIhX9RWw9fuhg621raxZO7gK3Hc4PMnT+MHHz+C1zbXcsldbw+5Qm+gNPAPUiC6nNOBGX8yi7f6sqt7BlPTb/fpyVQ5Z57PS2l+dr8z/te31mEMnDanHIBz51fwzvbI1o9q5Hltc6TH0xLr92E4u/K4Kfzus8ewu6GDi371Ju/tbkrbuRMt4PpEolvaRuhSTi/gOjDjTz7wTxqdz8JJpYNK99h9ejKV6gEYX5zL3n5m/Mu21DEqP5v5EyK1BeceMZ6wgRfXJ26O9e6OBs77+eus192+DimvbKpl6pj8mBdsh6Mlc8byxJdOJCfLw2W/fZtn1qSn1XOiGf8FCW5LUz80d+uOsYCre4ipnsKc/mv4+1o6v4K1e5rZWDOwFg6ZnvEDVJbmUd0UP/CHw4blW/Zz8qzyntK3ueOLmFZWwN/XxX+xW7O7ic/9/l021LRw5ysfOD5ulRldgRBvfVg/pAVXbjR7XBFP3XAS8yeUcOMjq9OyUDFugtcY87lUnlhEvMAKYI8xZli9kBhjIuWcVoonxxupQx/qjH/iqLyE7WBjWbqwgp/+YwsX/PINLlxUyfWnzmBOEiWaPdsuZmgBF8D4ktyEb2837m1hf1s3p80+8LZeRDj3iPH8dvk2Gtv9jCrovf3jxpoWPnPfO4wqyObcI8bz+KoqdtV3JF0mq9zr7Q/r6Q6GOX2Y5/djGVOYw4PXHseXHlrFrU+uwx8M87mTUle2nFSOX0TOF5Fvisht9s2Bc38F2OjA86RdIBRZeedsqidxH/54KkryePFrp3LVCVP4+9q9nPOz5Vxz/7vsa0mcQumZ8edncMZfcmARVyzLt0QWbZ06q6zX8fPmVxAKG17a0Dvd80FtG5++91/k+7w8fO3xfOOcOXg9wn1vbk/ND6DS6tXNteRle/utlx+ucrO93PXpoznn8HH85zMb+G0Ku3r2G/hF5C7gMuBGQIBLgSlDOamITATOB+4dyvNkil2941TgN8awJ8lVu7FMHJXP9y44nLe+fTo3nzWbZVvquP+tHQm/pyfHn6GqHoDxJQcWccWybEsth1UUM7ZPzf7hlcVMHJXH02uqWb2rkcdXVnH7C5v41D3/RER46NrjmDQ6n3HFuVywsJK/rNjd05FRDU/GGF7ZVMtJM8v6Xek9nPmyPNz5qaNYuqCCH/59E798eWtKzpPMjP9EY8xngEZjzH8CJwCzh3jenwHfBOJGShG5TkRWiMiKujp3tTS1A7xdzeP1CF6P4A8NrgFTS2eQ1u5gwi0XkzGqwMeNZ8xiyuh8dtUnLvFs7gzg9ciAepA7rTJqEVdfkX79jb3SPDYR4bz5FbzxwX4+/uu3uPnRNdy1bBujC3w8dO1xvS78XXvydDr8IR5+Z1fqfhCVch/WtVHV2MmSucO/mqc/2V4PP7tsEZ84cgI/fmkLK3c2On6OZP7q7Xq7DhGpBOqBisGeUESWArXGmJUi8pF4jzPG3A3cDbB48eLMdzWK0hP4sw7MPOx9dwdj9yBKOROZNDqfnQ3tCR8TadCWNeBrCk4an2AR19sf1hMMG06dXXbQ1wCuP3U6E0rzqCzNY3p5AZNH58dc/j6vslrYu+AAAB4RSURBVJiTZo7h/re28/mTp2VklaQaulc22WWch15+P5Ysr4fbL13IhYsqOXrKKMefP5m/gr+JSClwO7AK2AE8MoRzngRcKCI7gD8Bp4vIg0N4vrQ7EPgP/PP5sjyDTvXYHScHm+rpy+5TkmjLt0y2a7AlWsS1fEsd+T4vi6fEzueOKczh6hOncta8ccwoL0zY8+TaU6azr6WbZ9emp1ROOe+VTbXMHV9E5RDfFQ8nXo+krIIpmc3W/8sY02SMeZxIbn+uMea7gz2hMeYWY8xEY8xU4HLgFWPMpwf7fJlgp3T6Bv7uQQb+wdTwJzJ5dD6tXcGeC7ixZLIXvy3RIq5Vuxo5avIoR2boH5ldzqyxhdyzfLsr9z9ViW2obmHFjsZhv1rXTRIt4Drd+hi9aOt84IyRvoCru0+O374/2Bl/VWMHBVYQdMLk0ZF3DjsT5PlbugIZLeW0VZTkHXRxtysQYvPeVhZMLInzXQMjIlx7yjQ21LTw9b+s4Y9v72DlzkZH9i1WqdXU4ef6B1dQVpjDNSksbxxpEuX4TwNeIbJgqy8DPDHUkxtjXgNeG+rzpNuBVM+B/LgvyzPoXj1VjZ1MGEQNfzxTxkS6Xe5s6GDhpNKYj2nuDFBZkvm3zRUluQct4tpY00IwbFgwMfbYB+OiRRN484N6Xttcy5Or9wDgEbj57DncsGSmY+dRA9fY7ueJ1Xu4cGEl5UU5PcdDYcONj6xmX3M3f77++F5fU0OTaAHX96yPKV3INRwdqOrpfXHXHxxcVU+yG7Aka9LoSEBP1LytpTOYkf12+6qIsYjr/apIm4WFk5yZ8UOkRvoXVxyJMYaa5i7WV7fwlxW7uf2FzcwcW8g5h4937FwqeRtrWrjujyvY3dDJL17eyi3nzuWTiyfh8Qg/fnEzr2/dzw8/MZ8jJzt/gXMki/uXLyJfT/SNxpifOD+c4aHvAi77/qAv7jZ2cMxU536x831ZlBflsLM+dmWPMYYWF1zchUjg77sT15qqJsqLchifgp77IkKlVQ10yqwyLvvt23z9z+/x1JdPcnxjDpXYc2truPkvayjOy+LOTx3JH97eybefWMsTq/Zw9uHj+PVrH3LFsZO44tjJmR7qISfRlbMi67YY+CIwwbp9ATgq9UNzr3gXdweT6mnuDNDSFXTswq5tyuh8dsWZ8XcHw/hDYdfk+KH3Iq73q5pZMKEk5aWmudle7rrqaPJ8Xq77w8qeRW0qtcJhwx0vbOZLD61ibkURz3z5ZJYuqORP/3Y8/3vxfDbva+W/n93Iokml/MeFh2d6uIekRJut/6e1YGsicJQx5mZjzM3A0cCIfgnuu4DLvj+YGf+eQfThT8bkBIu43NCgzWZvyFJtVfa0dQf5sK7N0fx+4vPn8atPHcWuhg6+9qf3XLER9qHuvje3c+erH3DZ4kn86brje1ZmezzCZcdM5uWbT+Pfz5nD3VcdTU7WobtKN5OSqZUbB/ijPvdbx0asbgfr+O0+/ENdtdvX5DH51LR00R3juoMd+F2R6intPeNft6cZY2CBg/n9/hw3fQzfXTqPlzfV8kvt5plSrV0BfvXqB5wyq4wfXTw/ZmAvK8zhhiUzD2rVoZyTzNW9PwDviMiT1ucfAx5I3ZDczw7wOQ7U8Ttdw2+bMiYfYyLPP6NP7/IWF8347Ty+vYjr/arIhd4FE9IX+AE+c8IU3tvdxM9f3sKJM8dwzNRDsxFYpt3z+vbIRunnzM3oqvGRLpkFXD8ArgEardvnjDH/k+qBuVnfJm32/cHk+KsaO8nL9jK6T3vhobJr+WOle9yU6snzeRkVtYhrTVUzE0rzGFOY3tI9EeG/PnYEk0bn89U/vadN3VKgvq2b372+jfPmj2e+Q2s01OAktSzSGLOSSJuGJ4F6EdEcP71z/DleD4FBBP6d9e1MGZPv+Oxn8uiCnufvyw2dOaONL8mjxqrlX1vV7GgZ50AU5mTx88uPZF9LF9/561pd5euwX736IZ2BEF8/a06mhzLiJdOW+UIR2QpsB5ZZH/+e6oG5mZO9erbXtzPVWnDlpLJCH/k+L7saDm6HYM9m3TDjh0iXzprmLhrb/exq6Ejbhd1YFk0q5etnz+bZ92t4dGVVxsZxqNnT1MmD/9zJpUdPYubYQ2PbxOEsmRn/fwHHA1uMMdOAM4F/pnRULudU4A+Gwuxu6GBqmfOBX0QilT0xunS2dEVaFbjh4i5EunTubeni/T2RhVvpzu/3df2pMzhh+hj+4+n1bKtry+hYDhU//8cWEPjKmbMyPRRFcoE/YIypBzwi4jHGvEqktn/E8ofCiECWJ6plwyDKOaubugiEDFNTtC3g5NH5Mfv1NHcGyPd5E3a0TKfK0jwa2v28u70BgCMynP/1eoSfXrYIX5aHm/60OmZllEreB7VtPLayiquOnzKiumu6WTJ/+U0iUggsBx4SkZ8DiZu9H+L8wTA+r6dXXn4wF3e3W/n3VMz4warlbzi4PXNzpzsatNnsyp4XN+xlenmBK8Y2viSX/7t4Aev2tHDHC5szPZxhKxw23PrkWvJ9WXzpIzMyPRxlSSbwXwR0AF8Dngc+JHbjthGj2wr80bKtjVgGsgDIvvA6LUWBf8qYfLqDYWpbu3sdb3FBS+ZoFaWRwL9lXxsLM5jf7+vsw8dz1fFTuOf17by2uTbTwxmWfv/WDv61vYHbls5Le6WWii+Zcs52Y0zYGBM0xjwA3Al8NPVDcy9/KHxQn/iefXcHMOvfvr+dfJ+XsSnqOjjZumjct3VDZBMWd1T0wIG2DQDzM5zf7+vW8w9jzrgivvHoGur6vICqxD6obeP/nt/EGXPHcuniiZkejoqSqB9/sYjcIiJ3isjZEvFlYBvwyfQN0X38wYMDf84gAv+O/e1MGVOQsoUssfryt3YFWF/dkrJ3GYMR3YwtU6Wc8eRme/nlp46ktSvIzY+u0ZYOSQqGwtz86BryfF5++In5uljLZRLN+P8IzAHWAtcCrwKXAh8zxlyUhrG5ViDRjH8AF3h31HcwrSw1F3Yh0gbCI7Arqpb/LyuqaOsO8unjp6TsvANlL+LyeoR5Fe4K/ACzxxXx3aXzWL6ljnvf2Jbp4bjOe7ub+PO7u3rajwDctexD1uxu4r8/doS2XnChRO/3pxtj5gOIyL1ADTDZGHPwBqkjjD9Gjt/+PNnAb5dynntE6vrA+7I8VJbm9aR6QmHDA2/t4OgpozJaKx9LRUke44oNeT53NuW68rjJvLF1P//7/GZmjS3SbQAtgVCYLz24kmqr5cb0sgKOmz6ax1ZWccHCSpYuqMzwCFUsiWb8PWvWjTEhoEqDfkSsVM9AZ/xVjZ0EwyZlFT22yaPz2WkF/lc21bKroYPPnTQ1peccjFvOm8ttF8zL9DDiEhF+/MmFHFZRxA0Pr2KdteZgpHtubQ3VzV18/6LD+e7SeUwek89fV1czpiCH/7pIWyq7VaIZ/0IRabHuC5BnfS6AMcYUp3x0LuXExd3tKa7osU0Zk8+L6/cB8Ps3t1NZkstHXbjb1CmzyjM9hH4V5GRx39XH8PFfv8Xn7n+XJ790ouPttIcTYwx3L9/GjPICPn3cFDwe4fMnT6M7GCIUNuT73FNAoHpL1I/fa4wptm5FxpisqPsjNuhD7HLOgaZ6duyPBP4pKVq8ZZs8uoD6dj8rdzbw1of1XHXCVLJcsnBrOBpbnMvvP3cMXYEQn/39uzR3BGjrDvLmB/v55ctb+a+/baArMDIWfL29rZ711S1ce8p0PFGLGXOyvBr0XU7/dwbBHwxT1KfB2UBn/DvrOyjweSlPcW2zXdnz/b9tJDfbwxXHTkrp+UaC2eOK+O1VR3P1fe/wkTtepbkzQHSxT2O7nx9/cuEhX8lyz/JtlBX6+PiREzI9FDVAGvgHwR8M9+rFDwPP8W/f387UstSVctrsdxRrdjdxxbGTKc13tv3zSHXijDJ+9amjeGxlFYdVFHPUlFEsmlTKA2/t4CcvbWFuRRHXnTr8V6p2B0Pc+PBqDq8s4cbTZ/bM7Lfsa+XVzXV87czZPXslq+FDA/8gxMrx5www8O+ob+eINCxWmjT6QCrpGhde1B3Ozj58PGf3uV5y4+kz2by3lR/+fdMhUf1zxwubeXHDPl7csI8tta38+NKF5GZ7uff1beRkebjqBPeUBavkabJ3EGKXc3p7vtafQChMVWMn01LQjrmvkrxsygpzOGVWGbPGFaX8fCOdiHD7pQs4bHwxNz2ymg9qW5P+3u5giOv/uIK/vV+dwhEm760P9nPvG9u58rjJ3HLuXJ5bW8Nld/+TDdUt/HV1NZccPdHxDYRUemjgH4SE5ZxJ5Ph3N3QQSkMpp+3+zx3Djz+5MC3nUpDvy+Keqxfjy/Jw7QMr2NN08J4IsdyzfBsvrN/HD57dOKi9HZzU3BHg5kfXMG1MAbeefxjXnzaDuz59NFv2tnLBnW8QCIf5/MnTMjpGNXga+AfBHwof1NI42xvJfSbzB7ujp5QzPaWAR0woYWyRrp5Mpwmlefz2qqOpa+3m3J8t5/l1exM+fmd9O7985QOmlxdQ09zFM2syO+v/7lPrqGvt5qeXLeqp0Dnn8PE8+oUTKC/M4YIFlUwv1w1VhisN/IMw1AVc2/dHFlSlYuct5R6Lp47m2ZtOYcqYAr7w4Epue2pdzFJPYwy3PbWeLI/w0LXHMWdcEXcv35aWrR/DYcOHdW1sq2tjb3MXLV0BnlxdxdNrqvnKGbNYOKn3Cu8jJpTwxreW8BN9Bzms6cXdQUgU+LuTSPXs2N9OUU6W5kdHgKllBTz+xRO5/YVN3PP6dt7d0cjtlyzodWH/ubV7WbaljtuWzqOiJI/rTp3OzY+u4bUtdSyZ4/zF4bbuIG9sreOVTbW8urkuZtfRoyaX8sU4/fN1Hcjwp4F/gIwx+ENhcvr88ucM4OLujvr0lHIqd/Blebj1/HmcOKOMf39sDRfc+QaXHzOZb5w9G1+Wh+//bT2HVxbzGatC5oKFldzx4mZ+u+xDRwP/h3Vt3L1sG0+u3oM/FKYoJ4tTZ5dzyqwycrO9tPuDdPpDBEKGi4+eoAH+EJb2wC8ik4A/AOMAA9xtjPl5uscxWIFQ5O33UFI9O+rbWTRplPODU662ZO5YXvnGR/jFP7Zy/1s7ePb9auZVFlPb2s1vr1rcE2h9WR6uOWkaP3huI2t2Nx2Ubhmo1bsauWvZh7y4YR8+r4dLF09k6YJKFk8d5ZrtN1V6ZeJ/PQjcbIyZR2QT9xtExL3dufqwq3YGG/j9wTB7GjuZluJWDcqdinOz+X9L5/H8V09h4aRS/rmtgSuPm8yiPsH98mMnUZSbxd3Le7eBbu4I8EFtG3uaOqlv66bDH4x7LcAYww+f28jHf/0W/9zWwJeXzOTNb5/ODz4+nxNmjNGgP4KlfcZvjKkh0uIZY0yriGwEJgAb0j2WwbADe986fq9H8HqEQD85/l0NHYRN6vbZVcPDzLFF/OGaY1lf3cKc8QevryjKzebTx0/ht8s+ZPmWOjbtbeEfG2tZubORUJ/NYCaPzueOSxdy7LTRPcfCYcN3n1rHQ//axZXHTeY75x1GQY5mdlVERn8TRGQqcCTwrxhfuw64DmDy5MlpHVciPYE/6+Bl6j5v/xuu283ZNPArEUm4evtzJ07ld69v5zP3vQPA3PFFfPG0GcwaV0h3IExnIES7P8if393NZXe/zfWnzuBrZ83CK8I3H3ufJ1bv4YsfmcE3z5mj15NULxkL/CJSCDwOfNUY09L368aYu4G7ARYvXuya/e4OBP6D3yb7sjz9pnp6avi1lFP1Y2xxLj+7fBH1bd0smTs2bgvoq0+Yyn8/u4G7ln3Isi11VJbk8vKmWv79nDncsGRmmkethoOMBH4RySYS9B8yxjyRiTEMlj8UqcOOF/i7+wn82/e3U5KXzSgt5VRJOG9+Rb+PKcjJ4oefWMCZh43jW4+/z8aaFm5bOo9rdGWtiiMTVT0C/A7YaIz5SbrPP1TdcXL89rH+Zvy7GjpS3oNfjUxnHDaOF792GrsaOg66WKxUtExc1j8JuAo4XUTes27nZWAcg2IH9r5tme1j/eX4q5s6qSzJS8nYlBpd4NOgr/qViaqeN4hs3zgs9Z/jj7/7kjGGmuYuTps9vFv1KqWGNy3kHSB7Rh+rBjq7n1RPc2eADn+IylJtmKaUyhwN/APU74w/QaqnuqkLgMpSTfUopTJHA/8AxVvAZR9LNOOvtvqya+BXSmWSBv4BiteywT6WKPDXNFuBv0RTPUqpzNHAP0CJqnr6q+Pf09RFtlcoK8xJ2fiUUqo/GvgHqL8Zf6JePdVNnYwvycXjGbZFTUqpQ4AG/gFKlOPP6adXT02z1vArpTJPA/8ADaVXT3VTl17YVUplnAb+ARps4A+FDXtburSGXymVcRr4B8gfCiMCWTHy9InKOWtbuwiFjc74lVIZp4F/gPzBMD6vJ2Z/80QLuHoWb2mOXymVYRr4B6g7GI6Z5gG7qscQDh+8fYAu3lJKuYUG/gHyh8Ixa/ghat/dGLN+e/FWheb4lVIZpoF/gPzBcNxNqu0Sz1iBv7qpi6KcLIpzs1M6PqWU6o8G/gHy95PqsR/T156mTp3tK6VcQQP/ANkXd2PpmfHHCPw1zZ2a31dKuYIG/gEKhAY3469u6qJCK3qUUi6ggX+A/EkE/r79eroCIRra/UzQVI9SygU08A9QdxKpnr4dOrWUUynlJhr4Byipi7t9Zvw1zZHFW5rqUUq5gQb+AfIHk6jj7zPj32PN+CfojF8p5QIa+AcoUY4/J07gr7HaNYwr0Q1YlFKZp4F/gBKXc3p7HhOtuqmT8qIccrK8KR+fUkr1RwP/AA0mx1/d3Kn77CqlXEMD/wAlU84Za8avFT1KKbfQwD9AkVRP7JRNrMBvjNHFW0opV9HAP0D+YJjsrNibpWd7I8e7o1I9zZ0BOgMh3XlLKeUaGvgHwBgTacsc5+JuToyLu1rKqZRyGw38A2BftB1Ijt8u5azQwK+UcgkN/AMQCEV21hpI4K9utts1aKpHKeUOGQn8IvJREdksIh+IyLczMYbBsAN6vDp+r0fweqRXk7bqpi6yvUJZgS7eUkq5Q9oDv4h4gV8B5wLzgCtEZF66xzEYPYE/wUIsn7f3huvVTZ1UlOTh8cS+IKyUUumWlYFzHgt8YIzZBiAifwIuAjY4faJbn1zLO9sbHHu+QD85fvtrf353N69uqgUiF3fnTyhxbAxKKTVUmQj8E4DdUZ9XAcf1fZCIXAdcBzB58uRBnaiyNI9Z4woH9b3xHDl5FCfNHBP3619eMpPVuxt7Pp81rpCLFk1wdAxKKTUUmQj8STHG3A3cDbB48WIzmOe4YclMR8eUjH87dXraz6mUUgORiYu7e4BJUZ9PtI4ppZRKg0wE/neBWSIyTUR8wOXA0xkYh1JKjUhpT/UYY4Ii8mXgBcAL3GeMWZ/ucSil1EiVkRy/MeY54LlMnFsppUY6XbmrlFIjjAZ+pZQaYTTwK6XUCKOBXymlRhgxZlBro9JKROqAnQP8tjJgfwqG4yQdozOGwxhheIxTx+gMt4xxijGmvO/BYRH4B0NEVhhjFmd6HInoGJ0xHMYIw2OcOkZnuH2MmupRSqkRRgO/UkqNMIdy4L870wNIgo7RGcNhjDA8xqljdIarx3jI5viVUkrFdijP+JVSSsWggV8ppUaYQy7wD4eN3EVkkoi8KiIbRGS9iHwl02OKR0S8IrJaRP6W6bHEIiKlIvKYiGwSkY0ickKmx9SXiHzN+n9eJyKPiEiuC8Z0n4jUisi6qGOjReQlEdlqfRyVyTFaY4o1ztut/+/3ReRJESl12xijvnaziBgRKcvE2OI5pAL/MNrIPQjcbIyZBxwP3ODScQJ8BdiY6UEk8HPgeWPMXGAhLhuriEwAbgIWG2OOINKK/PLMjgqA+4GP9jn2beBlY8ws4GXr80y7n4PH+RJwhDFmAbAFuCXdg+rjfg4eIyIyCTgb2JXuAfXnkAr8RG3kbozxA/ZG7q5ijKkxxqyy7rcSCVau25hXRCYC5wP3ZnossYhICXAq8DsAY4zfGNOU2VHFlAXkiUgWkA9UZ3g8GGOWAw19Dl8EPGDdfwD4WFoHFUOscRpjXjTGBK1P/0lkF7+MifNvCfBT4JuA6ypoDrXAH2sjd9cF1GgiMhU4EvhXZkcS08+I/OKGMz2QOKYBdcDvrXTUvSJSkOlBRTPG7AHuIDLrqwGajTEvZnZUcY0zxtRY9/cC4zI5mCRdA/w904PoS0QuAvYYY9ZkeiyxHGqBf1gRkULgceCrxpiWTI8nmogsBWqNMSszPZYEsoCjgN8YY44E2nFHeqKHlSe/iMiLVCVQICKfzuyo+mcidd6um6lGE5FbiaRNH8r0WKKJSD7wHeC2TI8lnkMt8A+bjdxFJJtI0H/IGPNEpscTw0nAhSKyg0jK7HQReTCzQzpIFVBljLHfLT1G5IXATc4Ethtj6owxAeAJ4MQMjymefSJSAWB9rM3weOISkc8CS4ErjfsWI80g8kK/xvr7mQisEpHxGR1VlEMt8A+LjdxFRIjkpTcaY36S6fHEYoy5xRgz0Rgzlci/4yvGGFfNVI0xe4HdIjLHOnQGsCGDQ4plF3C8iORb/+9n4LIL0FGeBq627l8NPJXBscQlIh8lkoK80BjTkenx9GWMWWuMGWuMmWr9/VQBR1m/r65wSAV+64KPvZH7RuAvLt3I/STgKiKz6Pes23mZHtQwdSPwkIi8DywC/ifD4+nFejfyGLAKWEvkby7jy/lF5BHgbWCOiFSJyOeBHwFnichWIu9UfpTJMULccd4JFAEvWX87d7lwjK6mLRuUUmqEOaRm/EoppfqngV8ppUYYDfxKKTXCaOBXSqkRRgO/UkqNMBr4lWuISMgqz1snIs8MtOuiiPyHiHzDuv99ETnTgTHlicgyqwFgSlmdRr+UwudfKiLfT9Xzq+FDA79yk05jzCKri2UDcMNgn8gYc5sx5h8OjOka4AljTMiB5+pPKRAz8FsN3obqWeACq6WAGsE08Cu3ehurwZ6IFIrIyyKySkTWWg2wsL52q4hsEZE3gDlRx+8XkUus+zvsfugislhEXrPunxa1gG61iBTFGMeVWCtY441DRKZaewHcY/Xdf1FE8qyvHWP1jX/P6iO/zjp+uIi8Yx1/X0RmEVkwNSPqsR8RkddF5GmsFcki8nXrHdE6Eflq1Pk3WT/zFhF5SETOFJE3JdJb/1jo6b/zGpFWB2okM8boTW+uuAFt1kcv8CjwUevzLKDYul8GfAAIcDSR1bD5QLF1/BvW4+4HLrHu7wDKrPuLgdes+88AJ1n3C4GsPuPxAXujPo83jqlEmoUtsr72F+DT1v11wAnW/R8B66z7vyTSZ8Y+T571POuizvcRIo3nplmf2z9vgTXe9UQ6u9rnn09kMrcSuM8a20XAX6Oe80rgl5n+v9ZbZm8641dukici73GgJfBL1nEB/sdqy/APIu8ExgGnAE8aYzpMpLvpQPsyvQn8RERuAkrNgR7vtjIgur9/vHFApBHbe9b9lcBU6xpFkTHmbev4w1HP9TbwHRH5FjDFGNMZZ4zvGGO2W/dPJvLzthtj2og0fDsl6vxrjTFhIi8ILxtjDJEXiqlRz1dLpEuoGsE08Cs36TTGLAKmEAmydo7/SqAcONr6+j5gINsXBjnwu97zfcaYHwHXEpltvykic/uOp895Eo2jO+pxISLvDuIyxjwMXGid4zkROT3OQ9sTPU+U6POHoz4P9xlLrnVONYJp4FeuYyIdF28CbrYuapYQ2RsgICJLiLwwACwHPmZV3hQBF8R5yh1E0iQAF9sHRWSGNUv+XyKdXXsFfmNMI+CVA3vkxhtHvJ+jCWgVkeOsQz1bLorIdGCbMeYXRK4hLABaiTQfi+d16+fNl8iGMx+3jg3EbCLpJzWCaeBXrmSMWQ28D1xBZKONxSKyFvgMsMl6zCrgz8AaIrswvRvn6f4T+LmIrCAyG7d91bpI+j4QIPZOTi8SSbEQbxz9+Dxwj5XCKgCareOfBNZZx48A/mCMqSfyzmOdiNze94msn/d+4B0iO7bda/07DcQSItU9agTT7pxKJSAiRwFfM8ZcNcjvL7Ty8YjIt4EKY8xXnBzjAMYyDnjYGHNGJs6v3MOJ2mClDlnGmFUi8qqIeM3gavnPF5FbiPyt7QQ+6+gAB2YycHMGz69cQmf8Sik1wmiOXymlRhgN/EopNcJo4FdKqRFGA79SSo0wGviVUmqE+f8TiosCb/SoKAAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "plt.plot(irdf3.bins, irdf3.rdf)\n", "plt.xlabel('Radius (angstrom)')\n", "plt.ylabel('Radial distribution')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "If you are splitting a residue over your two selections, you can discount pairs from the same residue by choosing appropriately sized exclusion blocks." ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "THR has these atoms: N, H, CA, HA, CB, HB, OG1, HG1, CG2, HG21, HG22, HG23, C, O\n", "THR has 4 carbons\n", "THR has 10 non carbons\n" ] } ], "source": [ "first = thr.residues[0]\n", "print('THR has these atoms: ', ', '.join(first.atoms.names))\n", "thr_c1 = first.atoms.select_atoms('name C*')\n", "print('THR has {} carbons'.format(len(thr_c1)))\n", "thr_other1 = first.atoms.select_atoms('not name C*')\n", "print('THR has {} non carbons'.format(len(thr_other1)))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The `exclusion_block` here ensures that the RDF is only computed from threonine carbons to atoms in different threonine residues." ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 13, "metadata": {}, "output_type": "execute_result" } ], "source": [ "thr_c = thr.select_atoms('name C*')\n", "thr_other = thr.select_atoms('not name C*')\n", "\n", "irdf4 = rdf.InterRDF(thr_c, thr_other,\n", " exclusion_block=(4, 10))\n", "irdf4.run()" ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "Text(0, 0.5, 'Radial distribution')" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEHCAYAAAC0pdErAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3deXxcdbn48c8z2ZemaZt0T7pRKKXQhVAom7IKiIC4URFB0IJXL6D4Q/G6XNGreL2uIGIFLCAiKCCgyCJbWQrdNyh039IlSdckbZJZnt8f50wySc5MJu2s6fN+vfLKzDln5ny7ZJ58v9/n+3xFVTHGGGO68qW7AcYYYzKTBQhjjDGeLEAYY4zxZAHCGGOMJwsQxhhjPFmAMMYY4yk3WW8sIlXAg8AQQIHZqvprERkIPAqMBjYCn1bVPR6vvxr4jvv0R6r6QE/3rKio0NGjRyek/cYYcyRYtGhRg6pWep2TZK2DEJFhwDBVXSwi/YBFwGXANcBuVb1DRL4FDFDVb3Z57UBgIVCDE1wWASd6BZJINTU1unDhwsT/YYwxpo8SkUWqWuN1LmlDTKq6XVUXu48bgVXACOBSINwbeAAnaHT1EeBFVd3tBoUXgQuS1VZjjDHdpWQOQkRGA1OBd4AhqrrdPbUDZwiqqxHAlojnW91jxhhjUiTpAUJESoHHgZtVdX/kOXXGtw5rjEtEZonIQhFZWF9ffzhvZYwxJkJSA4SI5OEEh4dV9Qn38E53fiI8T1Hn8dJaoCri+Uj3WDeqOltVa1S1prLSc57FGGPMIUhagBARAe4DVqnqLyJOPQ1c7T6+GnjK4+XPA+eLyAARGQCc7x4zxhiTIsnsQZwGXAWcLSJL3a+LgDuA80RkDXCu+xwRqRGRewFUdTfwQ2CB+3W7e8wYY0yKJC3NNR0szdUYY3onLWmupoM/GOKxBVsIhfpOMDbG9H0WIFJg3rpd3Pr4cpZt3ZvuphhjTNwsQKRAiz/ofg+luSXGGBM/CxApEHCHlgIhCxDGmOxhASIF/EEnMASCNgdhjMkeFiBSwO8Ghrag9SCMMdnDAkQKBKwHYYzJQhYgUsBvcxDGmCxkASIF/AEnMLQFLEAYY7KHBYgUCPccArZQzhiTRSxApEB4kjpgk9TGmCxiASIFwmmubTZJbYzJIhYgUiBgPQhjTBayAJECfpuDMMZkIQsQKeAPuAvlLIvJGJNFLECkQEcWkwUIY0z2sACRAh1ZTDbEZIzJHrnJemMRuR+4GKhT1UnusUeBY9xLyoG9qjrF47UbgUYgCASi7XaULTqymKwHYYzJHkkLEMAc4C7gwfABVf1M+LGI/BzYF+P1Z6lqQ9Jal0JWi8kYk42SFiBUda6IjPY6JyICfBo4O1n3zyRWi8kYk43SNQdxBrBTVddEOa/ACyKySERmxXojEZklIgtFZGF9fX3CG5oI4R6E33oQxpgskq4AMRN4JMb501V1GnAh8BUROTPahao6W1VrVLWmsrIy0e1MiHBg8NschDEmi6Q8QIhILnA58Gi0a1S11v1eBzwJTE9N65LDdpQzxmSjdPQgzgXeV9WtXidFpERE+oUfA+cDK1PYvoQLWA/CGJOFkhYgROQRYB5wjIhsFZHr3FNX0GV4SUSGi8iz7tMhwBsisgyYD/xTVZ9LVjtTwd8+B2EBwhiTPZKZxTQzyvFrPI5tAy5yH68HJierXenQkcVkQ0zGmOxhK6lTIGA9CGNMFrIAkQJ+S3M1xmQhCxApYPtBGGOykQWIFAjvB2E9CGNMNrEAkQLh/SBsDsIYk00sQKRAwHaUM8ZkIQsQKWClNowx2cgCRArYQjljTDayAJECAdtRzhiThSxAJJmqRmQxWQ/CGJM9LEAkWTCkqNtxsDRXY0w2sQCRZOHMJRHbUc4Yk10sQCRZeFipOC8Hf1BRtV6EMSY7WIBIsvCwUlG+UzjX1kIYY7KFBYgkC9dfKsr3uc8tQBhjsoMFiCQL7wVRnJfrPrd5CGNMdkjmjnL3i0idiKyMOPbfIlIrIkvdr4uivPYCEflARNaKyLeS1cZU8AfCPYicTs+NMSbTJbMHMQe4wOP4L1V1ivv1bNeTIpID/Ba4EJgIzBSRiUlsZ1KFM5eK8nLc5zbEZIzJDkkLEKo6F9h9CC+dDqxV1fWq2gb8Bbg0oY1LofAkdXG4B2GL5YwxWSIdcxBfFZHl7hDUAI/zI4AtEc+3useykj/YZYjJJqmNMVki1QHid8A4YAqwHfj54b6hiMwSkYUisrC+vv5w3y7h2tNcw0NM1oMwxmSJlAYIVd2pqkFVDQF/wBlO6qoWqIp4PtI9Fu09Z6tqjarWVFZWJrbBCRAOCMXWgzDGZJmUBggRGRbx9OPASo/LFgDjRWSMiOQDVwBPp6J9ydB1oZzNQRhjskVust5YRB4BPgxUiMhW4PvAh0VkCqDARuB699rhwL2qepGqBkTkq8DzQA5wv6q+m6x2Jpu/WxaTBQhjTHZIWoBQ1Zkeh++Lcu024KKI588C3VJgs1GgvQfhdNZsiMkYky1sJXWSdWQx2RCTMSa7WIBIsshqrmC1mIwx2cMCRJJ1DDHZQjljTHaxAJFktlDOGJOtLEAkWUc1V8tiMsZklx6zmESkEvgSMDryelW9NnnN6jsC1oMwxmSpeNJcnwJeB/4NBJPbnL7H320ltfUgjDHZIZ4AUayq30x6S/qobluOWoAwxmSJeOYg/hFtYx/Ts0CXYn02xGSMyRbxBIibcIJEi4g0ul/7k92wvsIfDOETKMj1tT83xphs0OMQk6r2S0VD+ip/KERujo/cHAFsRzljTPaIqxaTiFwCnOk+fVVV/5G8JvUtgaCS5xPyfNaDMMZklx6HmETkDpxhpvfcr5tE5CfJblhf4Q+GyMv14fMJOT6xAGGMyRrx9CAuAqa4m/wgIg8AS4DbktmwvsIfVHLd3kOuT6wWkzEma8S7kro84nH/ZDSkrwoEQ+S58w95OT7LYjLGZI14ehA/AZaIyCuA4MxFfCuprepD/MEQeTlOHM7LsSEmY0z2iCeL6REReRU4yT30TVXdkdRW9SH+kLZnMOXm+KwWkzEma0QdYhKRCe73acAwYKv7Ndw9FpOI3C8idSKyMuLYz0TkfRFZLiJPikh5lNduFJEVIrJURBb29g+VSQLBUHsGU55PbIjJGJM1YvUgvg7MAn7ucU6Bs3t47znAXcCDEcdeBG5z953+Kc5Ed7QyHmepakMP98h4gaCSl+vOQeT6rNSGMSZrRA0QqjrLfXihqrZEnhORwp7eWFXnisjoLsdeiHj6NvDJuFuapdqCoU5ZTNaDMMZki3iymN6K81hvXQv8K8o5BV4QkUUiMivKNQCIyCwRWSgiC+vr6xPQrMQKBLVLFpP1IIwx2SFqD0JEhgIjgCIRmYqTwQRQBhQfzk1F5L+AAPBwlEtOV9VaERkMvCgi76vqXK8LVXU2MBugpqYm4349D4Qis5h8VmrDGJM1Ys1BfAS4BhgJ/CLieCPw7UO9oYhcA1wMnKOqnp+Wqlrrfq8TkSeB6YBngMh0bUGlKN8dYrI0V2NMFok1B/EA8ICIfEJVH0/EzUTkAuBW4EOqeiDKNSWAT1Ub3cfnA7cn4v7p4GQxuUNMPhtiMsZkj3gWyk0SkeO6HlTVmB/aIvII8GGgQkS2At/HyVoqwBk2AnhbVW8QkeHAvap6ETAEeNI9nwv8WVWfi/+PlFmcOQh3iClXaPVbgDDGZId4AkRTxONCnOGhVT29SFVnehy+L8q123BqPqGq64HJcbQrK/iDoY6Fcj4fTcFAmltkjDHxiWcldad1ECLyf8DzSWtRH+MPdS21YZPUxpjsEG+xvkjFOBPXJg5d01yt1IYxJlv02IMQkRU46xIAcoBKsnjSONWcIaZwFpNVczXGZI945iAujngcAHaqqg2kx8nv7igH4VpM1oMwxmSHeOYgNrnF+U7H6Um8gbNhkIlDINhloZz1IIwxWSKeLUe/BzwADAIqgDki8p1kN6yv8Ac1YojJehDGmOwRzxDTlcDkcME+d4/qpcCPktmwvsLJYrJaTMaY7BNPFtM2nPUPYQVAbXKa07cEQ4oqndJcrRaTMSZbxCrWdyfOnMM+4F0RedF9fh4wPzXNy27h3kLkjnLWgzDGZItYQ0zhndwWAU9GHH81aa3pY8LBoOuOcqqKW0rEGGMyVk/F+sxhCGcsRc5BgDP0FO5VGGNMpoo1xPSYqn66y0K5dqp6QlJb1gd0DDH5On13MpvS1ixjjIlLrCGmm9zvF8e4xsTgD3XtQYh7PEQRFiGMMZkt1hDTdhHJAeao6lkpbFOfEQjPQUQslHOOWyaTMSbzxUxzVdUgEBKR/ilqT5/SfYhJOh03xphMFu9+ECvcNNfm8EFVvTFpreojwoX5IneUc45bgDDGZL54Fso9AXwXZ0/oRe7XwpivcInI/SJSJyIrI44NFJEXRWSN+31AlNde7V6zRkSujud+maYji6ljR7nI48YYk8niCRDlqvpA5Bfg+aHuYQ5wQZdj3wJeUtXxwEvu805EZCDOFqUnA9OB70cLJJmsretCOetBGGOySDwBwuu392vieXNVnQvs7nL4Upzif7jfL/N46UeAF1V1t6ruAV6ke6DJeN0nqcNzENaDMMZkvljrIGYCnwXGiMjTEafK6P6h3xtDVHW7+3gHMMTjmhHAlojnW91jWSUQ6jLEFM5isl3ljDFZINYk9VvAdpwS35H7UjcCyxNxc1VVETmsX6dFZBYwC6C6ujoRzUqYbkNMOTbEZIzJHlGHmFR1k6q+CpwLvK6qr+EEjJHA4dSJ2CkiwwDc73Ue19QCVRHPRxKlgqyqzlbVGlWtqaysPIxmJV77JHVELSawISZjTHaIZw5iLlAoIiOAF4CrcCafD9XTdMxrXA085XHN88D5IjLAnZw+3z2WVdrnINzspbxcWyhnjMke8QQIUdUDwOXA3ar6KeC4eN5cRB4B5gHHiMhWEbkOuAM4T0TW4PRO7nCvrRGRewFUdTfwQ2CB+3W7eyyrtA8xuT2IXJ8tlDPGZI94FsqJiMzA2VnuOvdYXIWEVHVmlFPneFy7EPhixPP7gfvjuU+milbNNRkBoi0Qoi0YorQgnn9SY4zpWTw9iJuB24AnVfVdERkLvJLcZvUN4WylrqU2krGr3K9fWs0n7n4r4e9rjDly9fjrpjs5/VrE8/WAldmIQ1sKexCbdx9kbX0TwZCS40vNXhMra/dRlJ/DuMrSlNzPGJNasdZB/EpVbxaRZ/DeD+KSpLasDwh021GuYz+IRGtq8RMMKbuaWxncr7DnFyTAzY8upX9RHo9/+dSU3M8Yk1qxehAPud//LxUN6YvCcxAd6yDCtZgS34Nobg0CULc/NQEiFFI27zpAUJV9B/30L8pL+j2NMakVaz+IRe7316JdY2Jri7IfRDKGmBpbAwDUNbYAya/OvrOxpf3PN29dAxdMGpb0expjUivWEJPnVqNhtuVoz7pVc01iLaZmN0Ds3N+a8Pf2snnXgfbHc9dYgDCmL4o1xBTeavQr7vfwkNPniBE4TIdAKIQI7ZPGuUmsxdQU7kGkKkDsdgLE0UNKmbu6HlVFJDWT48aY1Oip1MYm4DxVvVVVV7hf38RZ2Wx60BYMtfceILk9iKYWtwfR2JLw9/ayZfcBfAJXnFTN1j0H2RjRozDG9A1xraQWkdMinpwa5+uOeIGgttdfguTtKNcaCLbPB6SyBzGsfxHnHDsYgNfX1KfkvsaY1Inng/464G4R2SgiG4G7gWuT2qo+IhAMtQ8rAfh8gk8SX4spnMEE4Unq5Nu8+wDVA4sZNaiE6oHFzF3dkJL7GmNSJ56FcouAySLS332+L+mt6iPagtppiAmcCetE9yDCE9R5OZLCHsRBzpng9B7OGF/B35fU0hYIkZ9rnUtj+oq4f5pVdZ8Fh94JBEPt8w5hToBIbA+i0Z1/GDWohPqmVoJJKOUR6UBbgIamVqoGFgFwxvhKmtuCLNm8J6n3Ncaklv26l0SBkLYvjgvLzZGEZzGFM5jGVpQQDCm7m9sS+v5dbd1zEICqgcUAnHrUIHJ8wutrbJjJmL7EAkQSdc1iguQOMY11ayLt3J/ceYjwGohqN0CUFeYxtaqcuTZRbUyfEmuh3OWxXqiqTyS+OX1LIBhqz1wKy/NJ4oeY2gNECQD1jcmdhwivgQgHCHCGmX710mp2N7cxsCQ/qfc3xqRGrEnqj8U4p4AFiB4Egl5DTL6E12IK9yDGpaoHsfsAJfk5nQLBGUdX8Mt/r+bNtQ18bPLwpN7fGJMasWoxfSGVDemL/CGvLCbBn+BJ5PAiuTEVTg+izqMH0RoIMn/Dbs4Yf/j7dm/ZfYCqgcWdVk5PHllOWWEur6+ptwBhTB8R1xyEiHxURG4Vke+Fvw71hiJyjIgsjfjaLyI3d7nmwyKyL+KaQ75fOvkDUbKYAontQYSHmMqL8hhYku/Zg3h8US1X3Te/Uw2lQxVeAxEpxydMHzOIRZssk8mYvqLHdRAicg9QDJwF3At8Eph/qDdU1Q+AKe575wC1wJMel76uqhd7HM8agVCofT/qMCeLKdEL5QKU5Ofg8wmD+xV49iBWbd8PwKbdzVQPKu52Pl6qyubdB/jQ0d17IuMGlzB3dX1KNy0yxiRPPD2IU1X188AeVf0BMAM4OkH3PwdY59Z86nP8QSUvN/lZTE0tAUoLnVg/uKyQOo8exOqdjQDUuimqh6q+sZXWQMgzyIweVEJbMMS2vYd3D2NMZognQIR/2g+IyHDADySqtvMVwCNRzs0QkWUi8i8ROS7aG4jILBFZKCIL6+szK83SHwx1qsUETj2mhAeItgClBW6AiNKDWFPXBHDYH97hDKaqgd4BAmCTFe4zpk+IJ0D8Q0TKgZ8Bi4GNRP9Qj5uI5AOXAH/1OL0YGKWqk4E7gb9Hex9Vna2qNapaU1l5+BOwieSdxSQJr8XU1NIRIIaUFVDf2EooYhiroam1ffHc1sMMEFv2uAFigEeAqHCObdjVfFj3MMZkhh4DhKr+UFX3qurjwChggqp+NwH3vhBYrKo7Pe65X1Wb3MfPAnkiUpGAe6aUPxRloVyis5haI4aY+hUSCCm7D3Sspg4PL+X65LCHmDbvcl4/ckBRt3ND+hVSmOdjU4MFCGP6glgL5c5W1Ze9FsyJSCIWys0kSk9ERIYCO1VVRWQ6TiDbdZj3Szm/50pqSXgWU3NrgEElzm/vQ8oKAGctREWp83jNTmd4qWb0AGoTMMQ0tKyQwrycbud8PmHUwBI2Wg/CmD4hVhbTh4CX8V4wd1gL5USkBDgPuD7i2A0AqnoPTqbUl0UkgDMHcoWqZt0udoGgkttlDiLX50t4LabGiEnqyn6FgLMvxHHucoTVOxspK8xlavUAFm7cc1hZRls8UlwjjRpUzHrrQRjTJ8RaKPd993vCF8ypajMwqMuxeyIe3wXclej7pppnFlOuL/FzEK0B+kXMQUDnfSHW7Gzi6CH9GFFeRCCk1DW2MKx/9yGieGzefYDTjoo+2jemooRXP7BUV2P6glhDTF+P9UJV/UXim9O3eGcxSfvub4mgqs46iIJwDyI8xNTafn51XSMXThrGCHfeYNveg4cUIFr8QXbsb+mhB+Gkum7fd5CRXSay2wIh6ptaGVF+aMHJGJNasSap+7lfNcCXgRHu1w3AtOQ3Lft13VEOEp/F1BoIEQhp+xBTQW4OA4rz2nsQ9U2t7D3g5+ghpe0fzFsPcaI6/LrqQdE/4MOZTF6prn98cwPn/PxVGlv8h3R/Y0xqRQ0QqvoDd2HcSGCaqt6iqrcAJwLVqWpgNvOuxZTYOYjwXhDhISZwMpnCPYjwBHV4iAk45InqLR5VXLsKr4Xwmqiet34XLf4QK2pt3yljskE86yCGAJE70LS5x0wP/FF2lGtLYBZTuFBfSWSAKOtYLBdOcR0/pJSSglzKi/MOebFcrEVyYUPLCinI9bGxy0S1qrJk814Alm2xAGFMNuixFhPwIDBfRML1ki4DHkhek/qGYEhRpXstJl9iazGFexClXXoQa+uc3d1W72yivDiPSjfldXj/okNeC7F59wEK83zt7+XF5xNGDSpmY5chpg0Nzew76AwtLd1iBf2MyQY9BghV/R8ReQ443T30BVVdktxmZb9wOY283C49iARnMXkFiMjV1Gt2NnL04H7tpblHDChi0yGuU9iy+wBVAzqX+fYyalBJt3uEew/HDiuzHoQxWSKuct+qughnUduTwC4RsTmIHrQHCI8d5dqCIRK1rCM8xBSepAanHlN4NfXqnY2MH1Lafm5EudODOJT7e5X59jJ6UDGbdh3oVO5jyZY99CvI5ZMnjmTH/hZ27EvupkbGmMPXY4AQkUtEZA2wAXjN/f6vZDcs24V7CV47yoEzBJUIzW1ePQhnsdyK2n3sbwkwfnBHgBg5oIjmtiD7DwZ6fa/aPQc9S2x0NbqihNZAiB0RVWUXb9rLlOpyplaXA7B0y95e398Yk1rx9CB+CJwCrFbVMcC5wNtJbVUf4HczlbyymICEzUM0tnjMQbiL5d5c48xDHD2kX/u54eFU1729q7h6oC1AY2uAIf0Le7y2aybTgbYA7+/Yz9SqciYOKyPXJyzbagHCmEwXT4Dwq+ouwCciPlV9BWdthInB7/YgumcxOc8TtViufQ6isPMkNcAba50AMT4iQIRTXbft7d0QT52bNht+71hGu1ufbmxwgtDyrfsIKUytHkBhXo47D2EBwphMF08W014RKQXmAg+LSB1gxXZ6EHADgFcWk3M+QUNMrQF8AkURxfPCq6nf39HIgOI8Kkrz28+FV1PX7uldDyKcNju4X/QMprBhZYXk5/raJ6rDE9RTqsrbvz+5pJZQSPFZOQ5jMlY8PYhLgQPA14DngHV4F/AzEdp7EB61mKAjgByuxhanzEZkZlFhXg7lxXmA03uIPDeoJJ+CXF+vF8uF97kOz2/E4vMJ1QOL24eYFm/ew5iKEgaUOIFqclU5Ta0B1tU39aoNxpjUimc/iGZVDalqQFUfwCmid0Hym5bdOrKYuu8oB4kbYmqOKNQXKfyb/tERGUzglGofUV7U+yGmXvQgwJmH2NhwoH2BXHhyGmBKVX8gMyaqQyFl74G2ni805ggUNUCISJmI3CYid4nI+eL4KrAe+HTqmpidOrKYutdiijx/uJoiCvVFCv+mHzlBHTZiQFGvd5ara2whP8fX3jPpyehBxWza3cyW3QdpaGplavWA9nNjK0rpV5CbERPVf3pnE6fd8bLVhzIJEQiGDnvPlUwSqwfxEHAMsAL4IvAK8CngMlW9NAVty2odWUzdS20ACavHFLmbXKTwPMT4wR4Borz3q6nr97dS2a+gx0VyYaMrSmjxh3j+3R0ATK3q6EH4fMIJVf0zYsHcv1fV0dwWZPnW9LfFZLcFG3fzsbve5PSfvsxzK3ekuzkJEStAjFXVa1T19zi7v00EPqKqS1PTtOwW3jXOa0c5gLZA4noQpTF7EKXdzg0vL6KhqZUWfzDu+9Q1tranz8YjnOr6xJJaCvN8TBjaOVBNHlnOqu37e9WGRGsLhFiwYTeQGcNdJjvV7W/ha48u5VP3zGPfgTaOGdKPrz+2lPd37E930w5brADR3udW1SCwVVUTtvxVRDaKyAoRWSoiCz3Oi4j8RkTWishyEcmqEuPhdQ5eO8o55xPUg2jxDhBXnFTF/3x8EoM86iaFU12392I18879LXHPP4CzsxzAqu37OWFkebehtslV5QRCyrvbov8Qrdq+n8vvfpMVSfrtftnWvRx0A5QFCHMo5q6u5+yfv8Y/l2/nq2cdxb9v+RAPXDud0oJcvvjAQnY3Z/f8VqwAMVlE9rtfjcAJ4ccikqjQeJaqTlFVr3UVFwLj3a9ZwO8SdM+U6KjF5J3F5E9gmqtXgBg1qIQrTx7l+ZqOVNf4h5nqGlvjWgMRNry8iHw3KEyLmH8IC6e8RlsPsXjzHj7z+3ks3ryXuWvq477vhobmuLOj3lq7CxE4e8Jglm7Zm7DyJ+bI0OIP8u0nVzC0fyHPf+1MvvGRYyjOz2VIWSGzP19DXWMr//HwovbPgmwUaz+IHFUtc7/6qWpuxOOyFLTtUuBBdbwNlIvIsBTcNyHa01w9ajE55xOU5hplkjqWjn0h4lsL0eIPsu+gv1c9iByfUDXQuU9kBlPYkLJChvUv9PzN/c21DXzu3ncYUJLPgOK8XqXDfvGBBZz/y7n86B/v0dwau5zIW+saOG54GWeOr6C+sbVXPSpj5ry1ka17DvKDS45jjLs4NGxKVTk//cTxvL1+N7c/817U9/AHQ/z5nc1s35eZE9txFetLEgVeEJFFIjLL4/wIYEvE863usazQvlAuSi2mRGQxhbcb7ecxSR3L0P6F+ARq40x1rXdTXONZAxEp/EMTOUEdafLI8m6ZTC+8u4Mv/HEBVQOK+esNMzh2WBnr6+Nbl7l930HW1TczfnAp976xgXN/8RrPrdzh2TM42BZkyea9nDqugiluD8dWd5t4NTS18tuX13LOhMFR92j/+NSRXH/mWB56exO3PLas2y8su5vb+Px98/n2kyu44U+LE7Y2KpF698mSWKeraq2IDAZeFJH3VXVub9/EDS6zAKqrM6fIrD8ULrXhPUntT8AcxEF/kJDiOcQUS16OjyFlhXEPMYW3L63sxSQ1wBnjK2nxhxgcJbBMrirnuXd3MOMnL3GgLcjBtiBtwRBTqsqZ84WTKC/OZ1xlKX9fWouq9phBNW/dLgB+/unJtPiD/NeTK7nhT4uYOb2an1x+fKdrF23aQ1swxIxxgzh2WD/yc3ws3bKXC4/Pmk6qSaNf/Xs1B/xBbrvo2JjX3XrBBArycrjz5TUs2bKHu2ZOY+LwMlZt38+XHlxIXWMrM6dX8cj8Ldzz2jq+evb4FP0J4pO2AKGqte73Onczouk45TzCaoGqiOcj3WNd32c2MBugpqYmYwaRO7KYvNNc/QnYVc5rN7l4DS8vinuIqaMOU+8CxNWnjubqU0dHPX/JlOGs2dmIzyeU5OdQlJ/LoJJ8Zp5c3R70xlaW0NgSoKGprT11N5q31u2ivDiPY4eW4fMJz0GnFO8AACAASURBVPzn6dz+zHs89PYmrjplFBOHd4yMzlvfQK5POGn0QApyczh2eJlNVPdBSzbvYXBZIcP7F8adot2T1Tsb+fM7m7nqlFEcNbh7lmCkHJ/w9fOOZsbYQdz0lyVcdvebfHZ6NY8t3EJpQS6PXT+DKVXlNLYE+PVLazhrwmCOG94/Ie1MhLQECBEpAXyq2ug+Ph+4vctlTwNfFZG/ACcD+1R1e4qbesjCWUpRF8oloJprY3g/6l4OMYEzDxHvB2LHKureDTHF04ZffGZKzGvGVTo/gOvqm2IGCFVl3rpdzBg7qL2+U16Oj2+cfwxPLqnld6+t486ZU9uvf2vdLiZXlbcHoikj+/PXRVsJhpQcqw/VJ7y9fhdXzHYKTw/uV8DU6nKmVg/g3GOH9PjBHsuPn11FSUEuN517dNyvmTFuEP+66Qxu+esy5ry1kSlV5fz+qhPbh21/eOkk3tmwm68/uoyn//M0CnJzenjH1EjXHMQQ4A0RWQbMB/6pqs+JyA0icoN7zbM4q7bXAn8A/iM9TT00HZPU3mmuPU1SN7cGetwzIjymWZJ/CAFiQBHb9x3stKlPNHWNLeT4hEEl+T1em2hjK515jJ7mIbbsPkjt3oOcOm5Qp+P9i/P43Cmj+OfybWxw98lubPGzfOu+TtdOqS7nQFuQNXWNCf4TmHT58zubKSvM5QeXHMdpR1XwwY5G7vjX+5z7i9e49Ldv8uC8jT2mofqDITY2NLNo026eW7mD37y0hlc/qOfGs8czsJc/D4NKC7j/6pN4dNYp/GXWKZ3m9AaU5PO/nziBD3Y28ssX1xzKHzcp0tKDUNX1wGSP4/dEPFbgK6lsVyK1p7l26UGEUz97SnP9yK/mcu6xQ/jvS46Leo3XbnLxGlFehD+o1DW2MrSHPR527m+lsrQgLZVXh/cvojDP12Mm01vrnNLmM7oECIDrTh/DH9/cwD2vruOnnzyBBRt3Ewxpp2snj3Q3Mtq8lwlDU5GkZ5JpT3Mbz63cwczpVc5Qp3t85/4Wnlm2jccX1/K9p97l9mfe48Ljh/GlM8ZwwsiOZIoWf5BHF2zh7lfXstMdYg2bOKyMz5/qnULeE59POHls9/+jAGdNGMwVJ1Uxe+46Jg4v46JJQ7uNQKRaOiep+7ToO8qFazFF70G0+INs3XOQh9/ZxKwzx7Zv8tOV137U8Qqvhdi650CPAaK3q6gTyecTxlSUsr7HALGLyn4F7UNSkSr7FfCZk6p4ZP5mbjp3PG+t3UV+rq/T+owxFSWUFTr1oa6YnjnJDubQPLGklrZgiJknd/63HFJWyBfPGMsXzxjLe9v287dFW3ls4RaeWbaNk8cM5ItnjGXb3oPtgWH66IHcct4xDC4roKI0/JWftA/u71w8kYWb9nDjI0v4Ub8CLp82kk/VjPT8f50KFiCSJNqOcrntWUzRexC73G6vP6jc89o6br90kud1hxMgwqUwNjQ0UzN6YMxr6/a3xLXVaLKMqyyJWStJVXlr3S5OO2pQ1InIWWeO5c/vbOYPr6/nnfW7qRnlbF4UJiJMripv37vCZC9V5ZH5m5lSVR6zNzhxeBnfGz6Rm88bz6Pzt3D/mxv40oNOUYfpowfyy09PYca46P+nkqG0IJd/3XQGr7xfx2MLt/KH19dzz2vrmFpdzuVTR3DxCcPby+anggWIJPEHvNNc8+PIYtrV5HRpR5QX8ZcFW/jqWUd5pop67SYXr6oBReTlCOsbel5jUN/YyrRR3VdDp8rYylKeXbGdFn+w04d62Lr6JhqaWpkRpesOMHJAMZdNHcGf39lMayDEN87vPsE4taqcu15Zy4G2AMWHMK9jMsOiTXtYW9fETz9xfM8XA2WFeXzpzLFcc9poXlpVR3lxHiePGZjSwBApL8fH+ccN5fzjhlLX2MKTi2t5YnEt333qXW7/x3ucdcxgLpg0lMlV5YwZVJLUoV/7KUiSQCiECN0yYnLjqOba4AaIWy84hq8/tozfz13Pdy+e2O26w+lB5Ob4GDWohHV1sYdu2gIhdjW39TrFNZHGVZYQUti06wDHDO1enfYtd/3DqeO8FyyF3fChcTy+eCsAMzyunVxVTkhhZe1+po+J3asymeuR+Vsoyc/h4hOG9+p1eTk+Lpg0NEmtOjSD+xVy/YfGMevMsby3fT9PLK7lqaXbeOG9nYDzs3/c8DImV5Vz24UTEh7U0jsD0of5g9qtzAZ0FO+LNUnd0OgMMU2rHsBlU0bw8Dub2oNGpKaWALk+oSD30P4Zx1aU9NiDCN830SmuvREef402D/HW2l2MKC9qL+0RzVGDS7lo0jDKCnM5YWT3XPPJ7orvpVv2HGaLTTIcbAv2WCF130E//1yxjUumjDik9UGZSkQ4bnh/vnvxRN6+7Wyeu/kM/veTJ/DxqSNoCYSYu7o+KT2evvM3mGH8wVC3RXIQsVAuxiR1Q7PzoVxRWsBXzhrHk0ucscjbLuy8arPZ3QviUP9jjK0s5ZUP6ggEQ1En3Xq7k1wyhFNdvTKZQiHl7Q27OO/YIXH9Pfz0kyewq6m129AfOH/fIwcUZcQ+Faa7m/6yhBdX7WT2VTWcN3GI5zVPLa2lxR9i5vQqz/N9QW6OjwlDy5gwtIxP1zh/znjS1Q+F9SCSJNqHbo5P8EnsWkwNjW3uyuIcxlaW8rHJw3lo3qZuOduNrYFDWgMRNq6yBH9Q2Rqj5EZdL/aiTpbi/FyG9y/0XAuxasd+9h7we6a3eiktyGXUoJKo56dUlduK6gz06gd1vPDeTkoLcrnxkSWsrO0exJ3J6S0cN7yM40dkzmrkVEjWPIQFiCTxh9SzBwHObwCxajE1NLV22sfhq2cdxUF/kDlvbuh0XVNL7wv1RRobsUo5mvYeRJrSXMPGVpZ6tjNcfyneANGTKVXl1O49mLHVNY9ErYEgP3jmPcZUlPDczWcysCSfa+csYFvE1p6tgSB3v7qOVdv3c8X06rRNMPc1FiCSxB8IeQ5jgJPJ5I+xo1xDUysVpR2pbOOH9OOUMYN4+YO6Ttc1t3nvBRGvcXGsUq7b34IIaVlFHWlcZQnr65u7VWZ9a90uxlaUMKx/YtJwzzl2CCLOKlyTGe5/YyMbGpr5/scmMqK8iPuvOYkDbUGunbOAptYAz7+7g/N/OZefPf8B5x47mE9My5qizxnPAkSSBELabZFcWG6OxMxi2tXURkWXneCmjSrn/e2NHGzr2KKzqaX3e0FEKi/OZ2BJPusbYvcgBpUUpH1F59jKUhpbA+2lxwH2Hmhj3rpdUcstH4oxFSWcP3EID87b5LmfRIs/yPefWsnGONKDzeHbsa+FO19ew7nHDuHDxwwG4Jih/fjtldNYU9fEh/73Fa5/aBH5OT4euHY69159kqUoJ5AFiCTxB0OeWUzg1GOKmcXU1EpFl0nhqVUDCISUlds6xl6b3EnqwzG2ooR1sXoQja0MSfPwEkQW7eto68PvbOagP8iVpyR25fP1HxrHvoN+Hl2wpdu5376ylgfmbeJ3r65L6D2Ntx8/u4pASPlelzTvDx1dyf9cNonCvBx+cMlx/OumM/jQ0ZVpamXfZQEiSZwspmhDTBI1iykQDLH7QBsVXYZ0plR31AoKa2oNUHqYvy2Nq4xdxqKusXd7USdL10ym1kCQOW9t5MyjKxNeO2la9QBOGj2A+97Y0OnfaUNDM79/bT35OT7+sXxbjzvWmcMzb90unl62jRvOHEu1u8d5pCumV/Pmt87m6lNHp72H21fZ32qSBIKxhph8UWsx7T7QhirdehAVpQVUDSxiSUSOflNLAnoQlSU0NLWx76Df8/zO/b3bizpZhpYVUpyf0z5f8tTSbdQ3tvKlM8Yk5X7XnzmO2r0HeXaFU2FeVfneUyspyPXxm5lTaG4L8s/lWVN9Puu88kEdsx5cyMgBRXz5w0eluzlHLAsQSeIPadTfanJzJGotpl1NTipr1zkIgClVA9prBYVCSnNb8LAmqaEjk8mrFxEMKbua0leoL5LPJ4ytLGFdfROqyh/mrmfC0H6cnsD5h0hnTxjMUYNLuee19agqz67YwetrGrjl/KP5yHFDGVdZwqMLuw9BmegWbdrDzv2xt7kN/9teN2cBVQOLefT6GRTlZ8beCEciCxBJEgiGyI/Sg8iP0YMIr1z2yhqaWlXO9n0t7NjXQnPboZfZiBRrv4VdTa2ElKhbhqba2IpS1jc08erqetbUNTHrzLFJS2f0+YRZZ4xl1fb9PP/uDn74j/eYOKyMz50yChHhMydVuTV/bP+IeGzdc4BP/34eF9/5Bu9u816I2BoI8v/+tpz/eXYVHzluKH/78gxGRKlkbFLDAkSS+IOh9s2BusrNkaiT1OEA0XWICWBqdUcpiOZWJ5vpcIeYqgcWk+sTzzUGmbCKOtK4ylK27jnIb19ey5Cygl7X2umtS6cOZ3C/Am58ZCk79rfwo49Pau8VXj5tJLk+8ZzINt09OG8TADkifOb3b/PW2oZO519fU8+ld73J3xZt5cZzxvPbz06zbKQMYAEiSfyx5iB8vqiT1LGGmCYOLyM/x8eSLXtpanXmDA63B5GX46N6ULFnDyI8HJApAWJsZQmqsHDTHr5w2hjyD7EGVbwKcnO49vQxtAVDXHFSVaf9IypKCzjn2ME8sbiWthiVef3BEFffP5+L73yd51buSFpJhEzW3BrgkfmbuXDSUJ78yqkMLy/k6j/O5+ll21i1fT+fv38+V903n6bWALOvOpGvn3d0WjanMt2lPESLSBXwIM62owrMVtVfd7nmw8BTQHjp8BOq2nXP6owWCIXaS3t35QwxeX9Q1De1kp/jo8yjZ1CQm8PE4WUs2byXC45zqk4eboCAjqGbrsI9iHSW2YgUTnUtyc9hZoo29fn8jFG0+kNcc+robuc+c1IVz7+7k5ff38kFk4Z5vv6Of73Pa6vrGd6/kBv+tIhjh5Vx0znjOX/ikCPmQ/CJxVtpbAlw7eljGNa/iL9efypfenAhNz6yBBGn3PZ3PnosV80YlTF7MRtHOvpwAeAWVV0sIv2ARSLyoqq+1+W611X14jS0LyH8gdgL5aL91tnQ2Mag0vyoY+tTqsp5dMGW9qyjwx1iAmeV8tzV9QRD2qk8ed3+jqKBmWBMRQmFeT4+e3I1/YvyUnLP4vxcbjp3vOe5M8dXMrSskEcXbPEMEP9Yvo373tjANaeO5jsfPZanl23jzpfXcsOfFjF9zEDuvbqGssLU/DnSJRRS/vjmRiZXlbf3wPoX5/HgddP50T/foyQ/ly9/eBzlxeldqW+8pXyISVW3q+pi93EjsAroc2vj/aHoFVKdWkzR5yBifSBPrS7noD/I4k1OuuvhFOsLG1dZSlswxNY9Bzodr2tsYWBJftKHcuJVlJ/D8zefya0XTEh3UwDn3/GTJ47ktdX17NjXOTtnbV0jt/5tOdOqy/n2RceSm+Pj8mkjefFrZ/KTy49nyeY9XHXvO+w74J1eHMvO/S3c/JclWbGa+7XV9axvaOba00Z3Ol6Yl8OPLjue2y461oJDBkvrT76IjAamAu94nJ4hIstE5F8iclyM95glIgtFZGF9fX2SWtp7gaDGGGKSqFlMu5o712HqamqV81vY6+4k3+EU6wuLlsnkrIHIjN5D2KhBJVEXIKbDp2uqCCl86vdv8dPn3mdl7T6aWgNc/9AiivNzuPvKEzsF2NwcHzOnV/O7K09k1fZGrrzvbfZ0qdLbk+89tZK/L93GrX9bnlFzGo0t3YPd/W9uYEhZARcd7z0EZzJb2n7SRKQUeBy4WVW77gKyGBilqpOBO4G/R3sfVZ2tqjWqWlNZmTlL7Z0spt5PUjc0dq/DFKlqYBGDSvJZ5pakTsgcRJSqrvWNLRmT4pqpqgcVc8/npjF6UAmz567n4jvf4JQfv8SGhmZ+M3MqQ/t7//2dO3EIv//8iaze2cRn732nfZvZnjy3cgfPv7uT6WMGMn/j7rjXYmzdcyDqYshEeGzhFk74wQtcee/bLN7s9G5X72zk9TUNfH7G6IwK6iZ+aflXE5E8nODwsKo+0fW8qu5X1Sb38bNAnogkZ0VUkjhZTNHTXL0mqVWVXc2dS313JSJMcbfGBBKya9bAknzKi/O67S5X15h5PYhMdMGkYTx03cks+K9zuePy4zl5zEB+eNmkHrdAPeuYwdz7+RrW1zfx2T/0HCQaW/x8/+mVTBjaj4e/eDKnjB3Ij59d1b5nRzSLNu3hvF/M5WN3vtGpRHaiPLW0lm8+vpzjR/Tn/e2NXH73W1w3ZwH/+9z7FOT6+GyKEgpM4qU8QIgz+3ofsEpVfxHlmqHudYjIdJx27kpdKw+fk8UUfaGc134Q+w768Qc15hATdKyHyM/1JWx+YFxlaaf9qfcd9FNvAaJXBpbkc8X0au675iSuPHlUXK858+hK/njNSWzc1cyV974Tc7jpZ89/QF1jK3d84gTycnz8+OPH0xoI8YNnuuZ3dHh/x36unbOAin757Glu44rZb1ObwCDx3ModfP2xZZw0eiCPzprB3FvP4v995BgWbNzNv1fVcfm0EQxIc6l4c+jS0YM4DbgKOFtElrpfF4nIDSJyg3vNJ4GVIrIM+A1whXbdCCDD+QOxJqnFcz+IBncNRGUPH8pT3WyQRAwvhUXuT72+vomP3/0mIs4HmEmuU4+q4N6ra1jf0Mzn7vOeuF60aQ8Pvb2Jq2eMZoq7d/bYylJuPPso/rliO/92N7GPtHnXAa66bz6FeT7+/MVTeOiLJ7PnQBtXzJ7XY5BYWbuPxxdtjRmwXv2gjv98ZDHHj+jP/decRFF+DiUFuXzlrKN4/Ztn86PLJvGN84/p5d+GySQpT3NV1TeAmAngqnoXcFdqWpQc/pj7Qfg894NoX0XdQ1rpCSP7I5LgAFFZyl8XbeW5ldu59W/LyfEJf7ruZE4em5id2kxsZ4yv5PdXncj1Dy7iqvvf4aHrTqZ/UR6tgSBb9xzk20+sYGhZId/4SOcP3FlnjuOZZdv57lMrGVSaz/DyIipKC2hoauXK+97GHwzx1+tnUDWwmKqBxfzpupP53H3vcMXsecz5wvT2tSVh727bx6/+vYYX3YCTlyOcPWEwH586klOPGsQHOxpZunkvS7bs4aVVdYwf3I8HvjC92//F/kV5fO6U+HpRJnPZWvYkcWoxxdhRzmMOor0OUw9DTP0K8xg/uJScKKU8DkU4k+mGPy3m6CGl3Hf1SVQN7F5i2STPWccM5nefm8YNf1rEhb+aC8D2/S2E+85/+HxNtw/i/FwfP/nE8Xzqnnl8/O63AGff8wJ36PHPXzqF8UP6tV8/uaq8PUic8/PXqOxXwKThZUwa0Z+1dU38a+UO+hXm8vXzjub08RU8u3w7f1+6jeff7dxDGTmgiI+eMIzvfHQi/Yv79lqOI5kFiCQIhpSQEr0Wk897P4hYZTa6+tq5R9MSCPZ4XbwmDisj1yeceXQlv75iCv36+AKuTHXOsUP4/VUnct8bGxjcr5DqgcVUDyxm4vAyjh3mve/FtOoBvHLLh1m9s5Ht+1vYse8gu5ra+FTNyPbhqEiTq8r553+ewUvv72Rl7X7e3baPuWsaKMrL4cZzxnPd6WPaFyJOqx7Aty6cwOtrG1ixdR/HDitjSlV5j8Ogpm+wAJEE4Q//2PtBePcgfAID4lg4dGGC88qrBhbz1rfOpqK04IgpAZGpzp4whLMnDOnVa6oHFXtuqhPr+i+c1rGXRovf+WWjMK97qYvcHB9nHTOYs9wtP82RwwJEEgTcHNRYC+X8oRCq2qmkRkNTKwNL8juVu0glW/Nw5PIKDMbY6pUk8Ad67kGoOkNRkep7WCRnjDGpZAEiCcJrHGKluUJHTyPMKbNhAcIYkxksQCRBeH4h1kI5oNtEtVOozxYVGWMygwWIJGifpI6RxeRc17kH4ZT6th6EMSYzWIBIgvAHf6w5CKBTRdcDbQEO+oM2xGSMyRgWIJIgvEo61kI5oNOeEA2N4TUQNsRkjMkMFiCSIFxnqadJan/ErnL14TIbtgDJGJMhLEAkQUcWUw9DTBH1mNrrMJVYgDDGZAYLEEnQkcUUfaEcdJ6kbi+z0c+GmIwxmcECRBJ0ZDFF31Eu8jqIKNRnPQhjTIawAJEEHbWYepiDiOhBNDS1UlaYm7ANgIwx5nDZp1ES9DzE1D3NtaGp1SaojTEZJV17Ul8gIh+IyFoR+ZbH+QIRedQ9/46IjE59Kw9dPNVcnesiexBWh8kYk1nSsSd1DvBb4EJgIjBTRCZ2uew6YI+qHgX8Evhpalt5eMLrG/KiBgh3iKlLFlOlBQhjTAZJRw9iOrBWVderahvwF+DSLtdcCjzgPv4bcI5E1sXOcOGho7weh5giF8q19riTnDHGpFI69oMYAWyJeL4VODnaNaoaEJF9wCCgIRkN+tidb7RvmJII+w46m873NEn9vadW8r/PvQ/A/paADTEZYzJK1m8YJCKzgFkA1dXVh/Qe4ypLaPPYAvRwVJYWMCzKBjxjK0qZOb2afQfb2o9NGFbGRccPTWgbjDHmcKQjQNQCVRHPR7rHvK7ZKiK5QH9gl9ebqepsYDZATU1N93084/CrK6YeyssOWX6uj59cfnxK72mMMb2VjjmIBcB4ERkjIvnAFcDTXa55GrjaffxJ4GVVPaQPf2OMMYcm5T0Id07hq8DzQA5wv6q+KyK3AwtV9WngPuAhEVkL7MYJIsYYY1IoLXMQqvos8GyXY9+LeNwCfCrV7TLGGNPBVlIbY4zxZAHCGGOMJwsQxhhjPFmAMMYY48kChDHGGE/Sl5YXiEg9sKmXL6sgSSU8EsjamBjWxsTIhjZCdrQzE9o4SlUrvU70qQBxKERkoarWpLsdsVgbE8PamBjZ0EbIjnZmehttiMkYY4wnCxDGGGM8WYBwC/1lOGtjYlgbEyMb2gjZ0c6MbuMRPwdhjDHGm/UgjDHGeDpiA4SIXCAiH4jIWhH5Vrrb05WIVInIKyLynoi8KyI3pbtN0YhIjogsEZF/pLst0YhIuYj8TUTeF5FVIjIj3W3qSkS+5v5brxSRR0TEe8ep1LbpfhGpE5GVEccGisiLIrLG/T4gA9v4M/ffermIPCki5ZnWxohzt4iIikhFOtoWyxEZIEQkB/gtcCEwEZgpIhPT26puAsAtqjoROAX4Sga2MewmYFW6G9GDXwPPqeoEYDIZ1l4RGQHcCNSo6iScUviZUOZ+DnBBl2PfAl5S1fHAS+7zdJpD9za+CExS1ROA1cBtqW5UF3Po3kZEpAo4H9ic6gbF44gMEMB0YK2qrlfVNuAvwKVpblMnqrpdVRe7jxtxPtBGpLdV3YnISOCjwL3pbks0ItIfOBNnnxFUtU1V96a3VZ5ygSJ3F8ViYFua24OqzsXZkyXSpcAD7uMHgMtS2qguvNqoqi+oasB9+jbOzpVpE+XvEeCXwK1ARk4GH6kBYgSwJeL5VjLwwzdMREYDU4F30tsST7/C+Q+e2E29E2sMUA/80R0Ku1dEStLdqEiqWgv8H85vktuBfar6QnpbFdUQVd3uPt4BDElnY+JwLfCvdDeiKxG5FKhV1WXpbks0R2qAyBoiUgo8DtysqvvT3Z5IInIxUKeqi9Ldlh7kAtOA36nqVKCZ9A+LdOKO41+KE8yGAyUi8rn0tqpn7lbAGfnbL4CI/BfOcO3D6W5LJBEpBr4NfK+na9PpSA0QtUBVxPOR7rGMIiJ5OMHhYVV9It3t8XAacImIbMQZpjtbRP6U3iZ52gpsVdVwD+xvOAEjk5wLbFDVelX1A08Ap6a5TdHsFJFhAO73ujS3x5OIXANcDFyZgXvaj8P5ZWCZ+/MzElgsIkPT2qoujtQAsQAYLyJjRCQfZzLw6TS3qRMREZwx81Wq+ot0t8eLqt6mqiNVdTTO3+HLqppxv/Wq6g5gi4gc4x46B3gvjU3yshk4RUSK3X/7c8iwifQITwNXu4+vBp5KY1s8icgFOEOfl6jqgXS3pytVXaGqg1V1tPvzsxWY5v5fzRhHZIBwJ6++CjyP80P4mKq+m95WdXMacBXOb+VL3a+L0t2oLPafwMMishyYAvw4ze3pxO3d/A1YDKzA+dlM+ypbEXkEmAccIyJbReQ64A7gPBFZg9PzuSMD23gX0A940f3ZuScD25jxbCW1McYYT0dkD8IYY0zPLEAYY4zxZAHCGGOMJwsQxhhjPFmAMMYY48kChMk6IhJ0UxdXisgzva3UKSL/LSLfcB/fLiLnJqBNRSLymlsIMqncyrT/kcT3v1hEbk/W+5vsYQHCZKODqjrFrXq6G/jKob6Rqn5PVf+dgDZdCzyhqsEEvFdPygHPAOEW+jtc/wQ+5paDMEcwCxAm283DLbQoIqUi8pKILBaRFW4xNNxz/yUiq0XkDeCYiONzROST7uON4Zr8IlIjIq+6jz8UsVhxiYj082jHlbgriqO1Q0RGu3tR/MHd9+EFESlyz53k7l2w1N3LYKV7/DgRme8eXy4i43EWpo2LuPbDIvK6iDyNu0JcRL7u9rBWisjNEfd/3/0zrxaRh0XkXBF5U5y9HaZDe32lV3HKVJgjmaral31l1RfQ5H7PAf4KXOA+zwXK3McVwFpAgBNxVicXA2Xu8W+4180BPuk+3ghUuI9rgFfdx88Ap7mPS4HcLu3JB3ZEPI/WjtE4heOmuOceAz7nPl4JzHAf3wGsdB/fiVNLKHyfIvd9Vkbc78M4BQjHuM/Df94St73v4lQDDt//eJxfDhcB97ttuxT4e8R7Xgncme5/a/tK75f1IEw2KhKRpXSUmn7RPS7Aj91yGv/G6VkMAc4AnlTVA+pUxO1t3a03gV+IyI1AuXbsMxBWAUTuLxGtHeAU5FvqPl4EjHbnUPqp6jz3+J8j3mse8G0R+SYwSlUPRmnjfFXd4D4+HefP26yqTTiF/86IuP8KVQ3hBI6XVFVxAsroiPerz43qHQAAAdtJREFUw6kqa45gFiBMNjqoqlOAUTgfxuE5iCuBSuBE9/xOoDfbdgbo+Jlof52q3gF8Eee39zdFZELX9nS5T6x2tEZcF8TpbUSlqn8GLnHv8ayInB3l0uZY7xMh8v6hiOehLm0pdO9pjmAWIEzWUqdK543ALe7kbH+c/Sn8InIWTgABmAtc5mYa9QM+FuUtN+IMzwB8InxQRMa5v3X/FKcScKcAoap7gBzp2EM6Wjui/Tn2Ao0icrJ7qH2rUREZC6xX1d/gzHGcADTiFKKL5nX3z1sszsZIH3eP9cbROMNe5ghmAcJkNVVdAiwHZuJsClMjIiuAzwPvu9csBh4FluHsLLYgytv9APi1iCzE+e0+7GZ3snc54Md7d7IXcIZ2iNaOHlwH/MEdOisB9rnHPw2sdI9PAh5U1V04PZmVIvKzrm/k/nnnAPNxdiG81/176o2zcLKZzBHMqrkakwAiMg34mqpedYivL3XnCxCRbwHDVPWmRLaxF20ZAvxZVc9Jx/1N5khEzrQxRzxVXSwir4hIjh7aWoiPishtOD+Tm4BrEtrA3qkGbknj/U2GsB6EMcYYTzYHYYwxxpMFCGOMMZ4sQBhjjPFkAcIYY4wnCxDGGGM8WYAwxhjj6f8DrRZZQOJjr0sAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "plt.plot(irdf4.bins, irdf4.rdf)\n", "plt.xlabel('Radius (angstrom)')\n", "plt.ylabel('Radial distribution')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## References\n", "\n", "[1] Richard J. Gowers, Max Linke, Jonathan Barnoud, Tyler J. E. Reddy, Manuel N. Melo, Sean L. Seyler, Jan Domański, David L. Dotson, Sébastien Buchoux, Ian M. Kenney, and Oliver Beckstein.\n", "MDAnalysis: A Python Package for the Rapid Analysis of Molecular Dynamics Simulations.\n", "Proceedings of the 15th Python in Science Conference, pages 98–105, 2016.\n", "00152.\n", "URL: https://conference.scipy.org/proceedings/scipy2016/oliver_beckstein.html, doi:10.25080/Majora-629e541a-00e.\n", "\n", "[2] Naveen Michaud-Agrawal, Elizabeth J. Denning, Thomas B. Woolf, and Oliver Beckstein.\n", "MDAnalysis: A toolkit for the analysis of molecular dynamics simulations.\n", "Journal of Computational Chemistry, 32(10):2319–2327, July 2011.\n", "00778.\n", "URL: http://doi.wiley.com/10.1002/jcc.21787, doi:10.1002/jcc.21787." ] } ], "metadata": { "kernelspec": { "display_name": "Python (mda-user-guide)", "language": "python", "name": "mda-user-guide" }, "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.3" }, "toc": { "base_numbering": 1, "nav_menu": {}, "number_sections": false, "sideBar": true, "skip_h1_title": false, "title_cell": "Table of Contents", "title_sidebar": "Contents", "toc_cell": false, "toc_position": {}, "toc_section_display": true, "toc_window_display": false } }, "nbformat": 4, "nbformat_minor": 2 }