Skip to content

Instantly share code, notes, and snippets.

@joezuntz
Created July 12, 2016 16:25
Show Gist options
  • Save joezuntz/84e431be1db26cc32d66be2847b12c30 to your computer and use it in GitHub Desktop.
Save joezuntz/84e431be1db26cc32d66be2847b12c30 to your computer and use it in GitHub Desktop.
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# PSFEx Im3shape Example"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"I'm going to assume that you have: \n",
" - an image file that you have already located your sources in \n",
" - a PSFEx output file, called \"psf_example.psf\" \n",
" \n",
"This first line loads the pylab interactive environment for this notebook. In a normal script to run on lots of images you would leave it out:"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Populating the interactive namespace from numpy and matplotlib\n"
]
}
],
"source": [
"pylab inline"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Load the python modules we will need"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"import py3shape\n",
"import numpy as np\n",
"import galsim\n",
"import galsim.des\n",
"from astropy.io import fits"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Use the example image file that comes from im3shape. It doesn't really make sense to use this with the PSF we load below because that is not really the correct PSF for this image."
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"image_filename = \"./example/image.fits\"\n",
"full_image = fits.getdata(image_filename)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Cut out a single object from the image to use. Because the example image file is just a regular grid we just pull out the top corner. In a real image you would need to find the object you wanted with x and y coordinates."
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAPwAAAD9CAYAAACY9xrCAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztnX+QVeV9xp+vgvEH/giSggXTZUY6mDEWC4lMg8lNS1Kc\nUIMjY0jHGTaVNCTBhHZoNR3a3CTMSB20TIdWOsHMEsmEWowESTVK652E6UCiYRVUFBK2ESKLBlGJ\nMcH49o97r2w5z8Pu3t2zrnmfz8yO1+9+7/uec+552d3n/Z7nGyklGGPy4JQ3+wCMMUOHF7wxGeEF\nb0xGeMEbkxFe8MZkhBe8MRnR8oKPiFkRsTsi9kTEjYN5UMaYcohW9uEj4lQATwGYCeAAgB8C+HhK\n6cnBPTxjzGAyosX3vRfA3pRSFwBExHoAHwXwxoKPCFf0GPMmklKKE2OtLvjxAJ7p8f/7AVxeyGpv\nrPkdVeCyav11JxnteTHLeSQ2U+Su7Pk/1cYXgHkkd/1+McgEHp5MYqN6vD5QBcY35nuV5KrpKiLO\nPpWeud+pAh9pzLdXjLGSxNh5qHitx+tfVoEzGvMdEWNgTzHUPomnni6GWL2v8WIlgMWN12N57vQz\ni7FxYlx2DwDAPDLfiIk8tyLGYJ/31p/z3Dnn1//7ZBW4uFp/re4N9vkB/FzYbbutsNYBtL7g+/bT\ne0e1/t+DNeDZGnBBpcXpjDEn5cUa8FKt17RWF/wBABf2+P8Lwf6tav5U31H1YjemTM6t1L+aHPgS\nTWtVpX8YwKSIaIuI0wB8DMAmmT2uIr9VDkM839lDPN+kIZ5vxBDPh+m/3fONqQztfD1oSaUHgIi4\nEvW/NE4FcEdK6eYTvp8wi4zdRgZbrWZ5icTO6c9hAnNIjP3dBQAXifgsElsgckeRmJqPHRvw//9+\nbqL+dlZ/l7O/6+TfziSmtJItIs5+V3xN5Kq/T7eR2Jh+5M4QuQ+LOPusLhG5THsC+N/gHSJ3FYnt\n6se4ANe76LWPQRXtkFK6D8B9rb7fGDP0uNLOmIzwgjcmI7zgjcmIlkW7XgeOSFQtGTe1GJsiBmGC\nlBJ8VNHFwceKsQ2X8txlYgwmlHSI3CUk1qmusVC1po8sxipiiOXdPD6HFKwogWgaiSmhS9QmoUYK\nb2aIwputr/D4MlJMs1nMx8S8gyL3YXGNcDaJ/YintglFsIuoa+vFRWLHp+45VnQGcPGWiY9dXLTz\nT3hjMsIL3piM8II3JiO84I3JCC94YzKiXJV+MRmbqaus5BAA2klMPUrbJeJsB2CjyFVjV0hMlVoy\ndVWVmG4+xuMbiUqvjm25iDNeI0o6ACwkaroqM66JeBeJsZJkQJ/LQhLbIHLXk5jaqVGl8mwnQpXW\nqs+QKe+qHJg9yrxtO8+dVXzaHAAvTV9HYket0huTPV7wxmSEF7wxGeEFb0xGlCvazSBjs9JAJaqs\nIWrE7Ot4rvJ2Y6IdK0UEdBkn8xFTwtNuElMijjgVWj7MxgWA2SLOrnNbP3LVMa/6qvjGlcXQHFFi\nWlPHQfwPLhL+B0wcVddzjYjvIiW+a0h5LwBsFWMw8XCRyGXHvFT43+EMEWdl1yRWs2hnTPZ4wRuT\nEV7wxmSEF7wxGdGypx0AREQX6k6TvwFwLKX03sE4KGNMOQxIpY+IfQCmppQOk+9x11rmxtklJmBK\nsXJerYg4K5VkjqeALpdlZaZqDJarxlW7BayMk5lUAPq82fEp99z+dMvZy5yEAawiaroy0VDHwRRv\n1TWGldGuFIr33PN5nH1W6v5SxiusfFg5DN9PYitE7hJRdt1OFHl2zKvLU+l5TxtjzLBjoAs+AdgS\nEQ9HxCcH44CMMeUxoL/hAbwvpfRsRLwDwIMRsTul9P03vrunejxzdAU4vzLA6YwxlAM14Ge1XtMG\ntOBTSs82/vtcRNyDehvp4wt+UnUgwxtj+sr4Sv2rySOD21sOEXFmRJzdeH0WgA8D2NnqeMaY8hlI\nb7mJAO5p/O8IAN/o2V8uIhKmkLGZqql6y7E6cVXDPlfEWS80VRet+rQxRVjVtjPbYdWPTZlMMIMO\npWyr39G2EJV3ClF4Ab5boHYWVE89VieuVPouEW8jMWVIwcr01Wei4kxNbxe5yjaa1ekfFbls50PV\n3atdIHbM7B4QtfQD6S23D9pR3hgzDHGlnTEZ4QVvTEZ4wRuTEQPdhz85rJcZa9GljCdYOaNyD62K\nOBN9VNlofwwpFO0kpgQfJRIyMa8mclW/uHFEoFOmFkygaxe5yvGXlXeKdmxSkGIiLXNkBbiQ2l9F\niZTtvuOzP6Wpz73nnXwMJtCpVXWQGG7sFYYbm4XD8DjiMKw+V4J/whuTEV7wxmSEF7wxGeEFb0xG\neMEbkxHl2lS3kbHnkGRlI8yU2C6R+xpRQAHgJqKCsn5gANAhTAcWEcVblU+yHYfnO3juJe08ziy3\nq2I+VTbKdgBUCTPbyVBGEBuUrTIxmVClw0pNZyXIKpft1mwV9/JKYdmwlMTU56pKYGlfRHEvVsi9\nqHaMVF++VeQNc0md8QbbVBuTPV7wxmSEF7wxGeEFb0xGlFtay55ZZmWmyrmzRmJ7WRDAzAqPM6dQ\nVfKJgvluncrYYkwdM31efybPXSzGqBDx6TohPKlzYc+4K+GJlecqz4A5wgGW5atjU8/a96fHHRUr\nxTVSpc1MoFO9+pSnARPoZolyWVZSrFagEg8xvhhqJ2kb+Lv9E96YjPCCNyYjvOCNyQgveGMyotcF\nHxFfi4juiNjZIzY6Ih6MiKcj4oGIULKIMWYY0WtpbURcgbpm+PWU0rsbsVsAPJ9SuiUibgTw9pTS\nTSe8L2EMGZspsUvE5KwU9Mg+njtvIo/THmT38dzKlTxeI+WkFwm1mp1LjafK3nJM1FequYqzclmh\n3NJyZ2XO0SHiXaTkcwzbpoHeLWBltHMe4bnTpxZj28S9oXZfZpEx1GeiSo3ZjoNyUO4iMWW6onYn\n2M4HU/qPtFha2+gk88IJ4asArG28Xgt+yxhjhhmt/g0/NqXU3XjdDYBsVBtjhhsDLrxJKaWI4H8X\n/KJ6/PXICnBaZaDTGWMYx2rAa7Ve01pd8N0RMS6ldDAiLgBwiGadVW1xeGNMvxhZqX81+dXg9pbb\nBGB+4/V8aC9TY8wwotef8BHxTQAfADAmIp4B8A8AlgO4KyKuR117vJa+mamu7IF/1TeNbfbNEGq8\nUkZZf7MxQo1Xm4sriCK/5CWe23lOMaaMIKh5AnjNteqxxswy1NjKLIMpzepaqM9qNlHklbmDEO+p\nqcUooqQDwsRE1LDTpnUA7id18PeLMbrE0Gw3Y+kmnjviqmKM7VoB8pAxncTYdRbj9rrgU0ofF98S\nT4QYY4YrrrQzJiO84I3JCC94YzKiXNdakLErJLk/LqZMtAC0qQITzJSIpspUVfkjYxmJqWNWBhHr\nSaw/JbQAN1BQ8zEBTJXWKsdfVgqqRELVW47NqdxbmdGIMq9QJcXrmANvN4kBmP4uHt9G7vGZwoiD\nHZ8yUlElvmytsPt5mV1rjckeL3hjMsIL3piM8II3JiO84I3JiFJtqiekPYXY/omTiomqXLOdxJTC\nvkbsNowgiqkq7VQwUwtpcvBEMVYVCu86MQbreVYVuWxXAOCGDcp4gpmEqNJaZvsNcLMSpeirHYDV\nDxZjr32I5+56rBjbdinPpdbh4CYme4WxieRbxdDWa3gqU9jVzonafZnGcskCEveFf8IbkxFe8MZk\nhBe8MRnhBW9MRpRbWjuKjM3sLteJYziPCG7zxISqtHYbeW59GnlmHeDPzgNcZKqIXFYSuULUxVaF\nYsMEMFXe259n3FV/tIUkxp5NB7QAxlC989Qz/Kw8typymdy8RuSqPm2srFWdnxqbXTs1HxNHFark\nlgnI7EH11S6tNSZ7vOCNyQgveGMywgvemIxotbdcNSL2R8SOxpd6atkYM4xotbfcFwG8nFK67STv\nS1hIxl5NDAZuEo1rlv9zMTbqczxXlSKyXQGlxiuFlhlYdIhchuqdp8pU2T+fqghaOd8yswWl6LNr\np0xJVLksU7GVIYU6F/a5bBH94pYT92K1C7GlxuPtlWJMlRSr3nJMkVf31wISU2XGLBfgu1Sb2TEM\nbm85ABC2HsaY4cpA/oa/ISIejYg73C7amLcGrT4tdzuALzdefwXArQCuL2T9sHr89e9WgPGVFqcz\nxpyUl2rAy7Ve01pa8CmlN3rJRcQaAPfSxPdUWxneGNNfzqnUv5o8O4i95RoNJJtcDWCnyjXGDB/6\notK/0VsOdQ/fL6JeST4FdR/qfQA+1aNffPN9CfjX4oDrPl2MXSeUWFbQPZMYaABaPWY9tpR6vEvE\nWe8uZVzArI8rwlRhnXD+WEAkYVXbriywWXtPZbjRQWJKlVbXmdWJq2uk6tKZYq36060jZhkrhFmG\nUpgWkM9qsfis1PVnj0l0igs9jzwQoe45pdKznSR2PQ9ylb7V3nJf6+19xpjhhyvtjMkIL3hjMsIL\n3piMKLm33KPFb0whzqLM+ADg5Z3MHAIAVoo4E32kIQJxQgWAReSYV90njuPKYkw5y6pSXiYyMZMD\nQDvwMt1IzddBYkqcm3dMfINQG8njM8Q9106KN7eK3A0kV90bW9Qxk+NTRiPr/p3HJ3ysGFNCKjP+\nUPe+EpYXk/uuSu65qg0wjMkeL3hjMsIL3piM8II3JiO84I3JiHJV+sVkbFZKqBThXeT9s8Rj+NvE\nGKwX1xZSlgkA00RpJjM0UD3u9vajtFaZdjA1XfXfE6ryubOLF/WMt71Cc+ehqEB/D1fQ3B8dYBcU\nwBZSi9vGU2m5M8B3M5Qt+U19fD+gy3MZzGYcAJYrpZ8Uq04W9+hucj+vEbkLxHyTyc4Cuxc3W6U3\nJnu84I3JCC94YzLCC96YjChXtLuIjM0eyFWiCnN1VY6gNRFnpYvqufCKiLMyxzaRy/qpKSFIiVdV\nElPPkFdEnDy/fevNn6Gpe4nq8wim0tw/E+ZGD5EDGQMiYAK46+/n0zh9Bl+VQbNSY9bnDdCfFXv+\nfhXpRQgAK0Q/QvYZLlM39PhiaJ4S7cQQbE2we+OIRTtjsscL3piM8II3JiO84I3JiJMu+Ii4MCIe\niojHI2JXRHyuER8dEQ9GxNMR8YAbURjz1uCkKn1EjAMwLqXUGRGjADyCere2TwB4PqV0S0TcCODt\nKaWbTngv7y3H3D/Vw/6sL5xQfrFMlK9WSUw5sjKFHeDGCqr3GutvphR2VS7LlF9R0SrHYEKx2C2Y\nPf8/CrF7D13Lk2/l4a/+Y7HG9y8fv5MnM0ddgO+oqLJrdv3p/QL9eS/vYwzQPeBqj5D5+A6H3FFh\n3L9HfIOU3Ha8qxhrb0GlTykdTCl1Nl4fBfAk6nsLVwFY20hbC32pjTHDiD7/DR8RbQAuA7AdwNge\nPvTdAET7V2PMcKJPraYav87fDeDzKaWXI47/ppBSSnX/OoJ7yxkzNDxZA3bXek3rdcFHxEjUF/ud\nKaXmX1/dETEupXSw0XbqEH2ze8sZMzRcXKl/Nfl2C73lov6j/A4AT6SUevrCbgLQrI+cDy3DGGOG\nEb2p9DMAfA/AY6j3kQOALwD4AYC7ALwTQBeAa1NKR054b8J0Mjaz6mXKNsDr45WSrmD5nf2YD+C2\nw2oMZhtdEbnK3GFm0aji9CNcjn9192gaHz3tQCG26tTP0tw/jm8XYreLQ6sK4R2k5d+fXM7r7v97\nu9iWYfeBsp5muxY1kbtExJfeXYwtuIbnshp2gNfvq2c12DGv205Tfy+dS+P/+yWyXVBlxiZn9b+3\nXEppK/RvAcop3RgzTHGlnTEZ4QVvTEZ4wRuTEeUaYLAeYszQQPRHm3RvsTfdnov/gCerElNWmvmq\nMDmYJkwO2voxH3MQVXsYqgcZ07RUeahwzx05oXiOF57/DM29J4r2uZf+RMz3Lzz87RUfLsTm/2ot\nyQReXDSOD7KZxJSwuYrElKrEyrkBfp3ZuAAXmwH+uSg3YnafKxMUVbtaIzGmxG2wAYYx2eMFb0xG\neMEbkxFe8MZkhBe8MRlRrkp/HRl7KUlWyiiJn37kME19dTIvMaUq/SIxnyqXZX3MlKnFRtK3bqvo\nWceuBQDUSKnkeWfy3IoYo60Y+vQ/3UZTj6BoWPQZIcdvxNU0/j/4o0KsWzw1/QrOoPFf/qp4ji9u\nFIo+U+/Z5wRoUwu2K6B2VNQOwLz7irEZV/LcdhLbIMZVhhtM1We9CN1bzhjjBW9MRnjBG5MRXvDG\nZES5ot1MMjZ7JliJdkeJQ+clI3muMspmJbCq15sqiayQmHrGuo3EmKgC6GPeTWLLRG6Nh0eSHml/\nfv43aO7b8OtC7Lv4U5pbwUM0/p/4SCH23IXv5Ae3kof58/CiDLoiyqAZ0nGWxNi1B/Q9yp6TV2W4\nTH/sErlKJGT381Y2rkU7Y7LHC96YjPCCNyYjvOCNyYhWe8tVI2J/ROxofM0amsM1xgyE3nzpjwH4\nq5695SLiQdQdbG9LKfFazSasDJC5kDJTDADYSBT5h086YxGmYKq+cCo+l/Szq4pedsxsgZVwAtr8\n4Pl9xVjXRJ4rPsFjXUUVe+MoXhb74kIiHxdbxQEA1k4Qkje5zufu5Y3hXlwuymXbiqGRYofjWI0E\n1W7IbOFWspe4V7D+dgCwiOwYAXzXaIYYg+0CLSafNQDMFJ83u5/ZR9LF396ba+1BNKrRU0pHI6LZ\nWw4ACpK/MWZ400pvuW2N0A0R8WhE3OF20ca8NehPb7kNqPeWOxoRtwP4cuPbX0G9ifD1hTceqB5/\nfXYFOKcykGM1xih+XgMO13pN609vuXXN3nIppUM9vr8GAG8xMr7al0M1xgyU8yv1ryY/HsTeco0G\nkk2uBrCz5QM1xgwZrfSW+zsAH0dd004A9gH4VI9+8c33JuAXZFRi5FAVB8B6dKnfSZTSz2yEVa4Q\nj6lxQYfI7SRq7nRR/7+tm8eZcYSo5R69sNhDDgAOV8cXg+r8KiR2urgvdgutlvVYU1bL7HkKAFhP\nYlzo589DqFw2LgCwzWRVd791k/gGM7tQN+kPiqG5l/PUDXt4fAVp4reE7CJhzKD2liM2H8aY4Y4r\n7YzJCC94YzLCC96YjCjXAGM2GZuVP6pKfCbaqdxVrOYQoHWOSsRhBgwAN7VQYzChSplXLBBxdo7K\neVU57bIxVI81Vgq6S+Sq8lVm+qCMRtTY7STWJnL3k5hyEmb3EcCvqRIalakFK/VWJdrs81ZOu6Lf\nIhVemQi60AYYxmSPF7wxGeEFb0xGeMEbkxFe8MZkRLkq/QIy9ppqMVYlMQCoPkaC3+K5c8UYG5ic\nKyRQZX7QRmLTRS7rFaaUX2ZxDHD76tn9yAW4UYIqPWVKs9r0WCHibAfgIpGrxm4jMWV4wj6rLpEr\n/C+o4q3MK9QY5LN9xzU/panPrSW23duKIQC6/JsdBzNYedUqvTHZ4wVvTEZ4wRuTEV7wxmREuaLd\nJWRsVqKoRK31pCZyjbBTVQIRK11Uwgx7dh7o3zHPY0HheLpCPCffQWLV/swHXmaqRDs2nyrDZa7D\nAC9fVc/fq7JRJlSp+VhZrBIUlfA3t4/HcDLYPVMRuYvJfTBP3ANMawaArcRDYRTxTzhq0c6Y7PGC\nNyYjvOCNyQgveGMyojfX2tMjYntEdEbEExFxcyM+OiIejIinI+IBN6Iw5q1Bryp9RJyZUnolIkag\nXhS5BMBVAJ5PKd0SETcCeHtK6aYT3pfQQcZm7qtKHWeKqfqnpSbii0hMmWXMEvI9M3JQpa5LSUy5\ntO4WcVa2q8wkVDlw9REy7lSey3YchDHDKQeZEzHw+kVnFYNL+BiypJjtDKj+e4tJrCZyN4hdksVE\nIVeGIjUhmy8gWw5rWEk4gHmXFmPrlXMxc6IFMO1dxVgbydvQokqfUnql8fI0AKcCeAH1Bb+2EV8L\nvUlljBlG9LrgI+KUiOgE0A3goZTS4wDG9vCh7wY1UjfGDDd6bTWVUnodwJSIOBfAdyPigyd8P9Wb\nThDuqR5/PbkCXFxp/UiNMZpDNeC5Wq9pfWomCQAppRcj4jsApgLojohxKaWDjbZTh+ibrq72dXhj\nzED4nUr9q8mTrfWWG9NU4CPiDAAfArADwCYA8xtp8wFsHODhGmOGgN5+wl8AYG1EnIL6Pw53ppT+\nKyJ2ALgrIq5H3XbgWvpuZmnMVGyljDI745kit9aPXmj3CzVeSY/MdOASkbuSxJgpBqDV+w4S669N\ndRtR5JXZAqspFzXsr68iajwAHH2pGJtxDs9Vx7yRKNO7zue57D5Qd7N6ZkFdU8ZW8QAAu59HEDUe\n4PfRDCF/TRFxdo4rVeF9397+BimlnQD+kMQPQy89Y8wwxZV2xmSEF7wxGeEFb0xGlGuA0UbGZuWT\nqjyUGTYoowtVplojMWWqoMpXmckEE2sA7pyqTBWUEQc7DmVeoYwq2DXtjwOsUnfUMTPdSH1WSqxk\nx6f2f9jYUsAUcVbSvUGsh+VE/AV4CbLqGcjKhFUfOib+qvzJJLbKBhjGZI8XvDEZ4QVvTEZ4wRuT\nEV7wxmREuSo97ix+o0JsppUizOLM0ALQKjbrb0d9pwHMm8TjrERYlday3JrIZeoqwM9b9WNb1o94\nm8hlJZ9KPV63ncdnXN73+ZSpBdtpEa7k9DqLj1XCVP2Dd/Pc2dfwOLvvVKUru6bK/KVLxNmOyhZm\n8HGaVXpjcscL3piM8II3JiO84I3JiHJFu3YyNus3pkpdmTOsEjmUEMTmUw/2qufFmcuqcq3dTGLq\nOXtVyttRDI3uPEBTD68bz8dg5bzquXwmCKrzY8+sA8Al5Ln1/l5nViY8S+TWSKw/PeQAYH0fjwEA\nFoo4K/3tjxCnrpG6Z9i1YyJvh0trjckeL3hjMsIL3piMaLXVVDUi9kfEjsaX+kvLGDOM6M3T7tWI\n+GDPVlMRMQNAAnBbSum2ITlKY8yg0JdGFKzVFAAIR4AedJD+ZjOJm6roY4ab9hVjcyfyXGV+wAwb\n5oidiTXilJjquuwVEgTQRRxS24VrqlJziVvv4elCjVefIOvBp34PY4pwTeRuEy6yzLCB9REEeF84\ngJf4KuV9KyknXSquM3M/BrgRR5vIZbs9AL/+0kT2iWJoL+kVBwCLyL0PAOvI/c92qDr421ttNQUA\nN0TEoxFxh7vHGvPWoC/NJF9PKU0BMAHA+yOiAuB2ABNRfxzgWQC3lnmQxpjBoZVWU9NSSrVmPCLW\nALiXv+vferyeCm1mZowZEE/XgD21XtNOuuAjYgyA11JKR3q0mvpSs69cI+1qADv5CJ/q+wEbY1rn\n9yv1ryb38d5yrbaa+npETEFdrd8Hr2xj3hKUW0vP1HBmXKDq4HfdR4If4LmzzuRxZmdcE/Opf/7Y\nDoCyEd5CYkrSZHX3AHCE7ABMEOenFG9mBqF2MpipiLLhVnXprMZb7ULURJyp+itELruP1Oenjpld\nD/V8wxahmi8gqrnaWegk/ffaRf899dxDX58NWe1aemOyxwvemIzwgjcmI7zgjcmIckW7yWTs3USQ\nWiIEKWbioAQwJdgw8UP1PKuIOBOOlGstc6JlQh6gzQ8Yq0Qp71xx7dh5V8TYNRJrE7lKaLykmwTF\nh7JUlOcyQwrVd5CVdCixTPUdZKXGSmhUPe6YYYYS/tpJTLnyKrGSjc1KeTdatDMme7zgjckIL3hj\nMsIL3piM8II3JiP6/LRcSzBDg9lEVVZlqu0k1iFymdEFwEswt7ADAzBXSLQHyW7DcmGW0U4k040T\neK5SYln6DKHGKxWb7UQotZqVtKqdBdXLbt7YYkzZUYvLgb1E6b+IjAsAq24vxhZ/mueqY55JTDQW\nCRONI6x/G4DFJF+dN9s5UYq+sstuIzHVd5Dgn/DGZIQXvDEZ4QVvTEZ4wRuTEeWW1q4mY9MeXTU+\nyORKMcae8wb0c+HMsVSVtC4RcZZf60euLNckrr4AMIc4+6qSVvYsO8DFIPU8/GZWtitEQuXeOp3E\nlPAkNFMqISvxkIm0Sryq9mM+JfCxMm+Ai6Oq7Jqdy9Gv8twpn+TxNhLbyNbwKS6tNSZ3hmbBP1Ub\nkmne4KdDPN/rv+XzyV9nSuLnQzzfj4d4vteGeL4eDM2Cf7o2JNO8wTNDPN9xE9/fzvmGesEfHuL5\nfjLE8/1miOfrgX+lNyYjvOCNyYhyVXpjzJsGU+lLW/DGmOGHf6U3JiO84I3JCC94YzKi9AUfEbMi\nYndE7ImIG4dgvq6IeCwidkTED0oY/2sR0R0RO3vERkfEgxHxdEQ8EBHKW3ew5qtGxP7GOe6ICOa/\n2spcF0bEQxHxeETsiojPNeKlnN9J5ivr/E6PiO0R0RkRT0TEzY14Ween5ivl/PpESqm0LwCnol79\n3gZgJOrV3BeXPOc+AKNLHP8KAJcB2NkjdguAv228vhHA8pLn+yKAvy7h3MYBmNJ4PQrAUwAuLuv8\nTjJfKefXmOfMxn9HoG5VMaPkz4/NV9r59fZV9k/49wLYm1LqSikdQ915/KMlzwkAwo5m4KSUvg/g\nhRPCVwFY23i9FsCckucDSjjHlNLBlFJn4/VRAE8CGI+Szu8k8wElfYYppeaTQqeh/gPpBZT7+bH5\ngBLv0ZNR9oIfD+CZHv+/H8c/0LJIALZExMMRIR45GnTGppSa/kzdAIQv06ByQ0Q8GhF3DOafEE0i\nog313yy2YwjOr8d8TYOoUs4vIk6JiE7Uz+OhlNLjKPH8xHxAyZ+fouwF/2Zs8r8vpXQZgCsBfDYi\nrhjKyVP997eyz/t2ABMBTAHwLIBbB3PwiBgF4G4An08pvdzze2WcX2O+DY35jqLE80spvZ5SmoK6\ns977I+KDJ3x/UM+PzFdByZ/fySh7wR8AcGGP/78QvDHOoJFSerbx3+cA3IP6nxVl0x0R4wAgIi4A\ncKjMyVJKh1IDAGswiOcYESNRX+x3ppSaDZZKO78e861rzlfm+TVJKb0I4DsApmIIPr8e800bivNT\nlL3gHwYwKSLaIuI0AB8DsKmsySLizIg4u/H6LAAfBrDz5O8aFDYBmN94PR+6E9mg0Lgpm1yNQTrH\niAgAdwD6WomkAAAAwElEQVR4IqXU00u4lPNT85V4fmOavz5HxBkAPgRgB8o7Pzpf8x+XBoN2fn2i\nbFUQ9V+tn0Jdrf9CyXNNRH0noBN1r5tBnw/ANwH8DMCvUdcnPgFgNOp+Jk8DeADAeSXO9xcAvg7g\nMQCPon5zjh2kuWYAeL1x/XY0vmaVdX5ivitLPL93A/hRY77HAPxNI17W+an5Sjm/vny5lt6YjHCl\nnTEZ4QVvTEZ4wRuTEV7wxmSEF7wxGeEFb0xGeMEbkxH/B4WNLwsljMK0AAAAAElFTkSuQmCC\n",
"text/plain": [
"<matplotlib.figure.Figure at 0x1128e3190>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"stamp = full_image[:40,:40]\n",
"pylab.imshow(stamp, interpolation=\"nearest\")\n",
"stamp_size = stamp.shape[0]"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Now read the im3shape parameter file."
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"options = py3shape.Options()\n",
"options.read(\"./example/example.ini\")\n",
"options.stamp_size = stamp_size"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The PSFEx files do not contain images of the PSF but instead decsribe the variation of the PSF across the image using basis shapes and polynomial coefficients. Here we use galsim to load those coefficents:"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"psf_filename = \"psfex_example.psf\"\n",
"psf_data = galsim.des.DES_PSFEx(psf_filename)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Use the position in the image and the size of the stamp to work out some parameters for the shape and location of the PSF. x and y here are the position of the postage stamp in the wider image. "
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"x,y = 20,20\n",
"nx = ny = (stamp_size+options.padding)*options.upsampling"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Get an array of the PSF at this position in the focal plane"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"<matplotlib.text.Text at 0x1133a2d50>"
]
},
"execution_count": 8,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAQMAAAEKCAYAAAAW3jADAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztvX+sZdd1HvZtzXskZ8gZjocUh+IPe2RTkkVLNWXLdFo7\nsV0YjoUEVlIEdoK2tmMjberGLtqitRwUtV0jhmPATgsXCVz/ktTYigUEdpUCSSQldmGksRUqoiWX\nokPaGog/xKFJejQzmSH57vD0j3PW3O9+71t77/t+zLtveBdwce89Z/9YZ++9vrXWt889twzDgLWs\nZS1recNBK7CWtaxlNWQNBmtZy1oArMFgLWtZyyRrMFjLWtYCYA0Ga1nLWiZZg8Fa1rIWAGswWMta\n1jLJGgwOqZRSzpZSLpdSLpZSniul/Eop5dbp3FeVUj5aSnmxlPKnpZRHSinvmc59cynltalevP6v\npI/3l1J+4npe11oOTtZgcHhlAPAXh2E4DuBrALwbwP80nfsnAP45gNMA7gLwQwAuUN1nhmE4Tq/3\nVvpY35X2OpE1GNwAMgzDswD+GYCvKqXcAeAMgF8YhmE2DMPWMAz/7zAM/2qHzRcAKKWcmSKK7y2l\nfH6KOv5mKeXrSimfniKQn7tWqZSvKKX8y1LKC6WUPyml/MNSyu10/mtKKZ8qpVwopXy4lPLrHIWU\nUv5iKeXRqd1/VUp55w71X0unrMHgcEsY6v0A3gPgU8MwvAjgSQC/Wkp5bynl9B73+TCABwD8VQD/\nG4C/DeA/BvBVAL6zlPLnqOzfAfAmAG8HcD+AH5v0vQnAbwD4ZQBfAuBDAP4SpiiklPIuAL8E4G8A\nOAXg5wF8ZKq3ln2SNRgcXikAfrOU8qcAfgfAbwP4yenctwA4C+BnADxbSvl/SikPUN17Jo8br7+y\nRL8/MQzDq8MwfAzARQC/NgzDC1N08jsA3gUAwzD80TAM/2KKTF4A8PcAfNPUxp8BcGQYhp8bhuHq\nMAy/AeAT1Md/AeDnh2H4N8MoHwTwylRvLfskGwetwFp2LAOA9w7D8C+3nRiGZwD8IACUUu4D8H8A\n+CCA/2gq8uwwDPfvsN9z9PmK+X7b1O9pjJHDNwI4jtHxvDSVuwfAM9LuU/T5ywB8dynlB+nYJsYo\nYy37JOvI4AaXYRieBvD3AbzjenU5vf8kgKsA3jEMw+0A/nPM19sXANwr9b6UPn8ewN8ZhuFL6HXb\nMAy/vp+Kv95lDQY3mJRSTpZSfnwi8N5QSrkTwPcB+Nc7aW4X5W8D8O8BXCil3Avgf6Bz/xrA1VLK\n3yqlbJRS3gvg6+j8LwD4m6WUh8sot5ZS/kIp5bYdXMNaOmUNBjeevIoxzP44gC8C+AzG8P17qUzv\ndqFuLbbq8fkfx7jl+UWMW53/OM4Pw/AqgP8EwPcD+FMA/ymA/3vSHcMwfBIjefi/Y0wtngDw3Z06\nr2WHUtYPN1nLKkgp5fcA/P1hGD5w0Lq8XmVfIoNSyreXUh4vpTxRSvnh/ehjLYdbSil/rpRy95Qm\nfA9GTuOfHbRer2fZ892EUsoRjOHdt2JkjP9NKeUjwzB8dq/7WsuhlrcB+DCAWwH8EYC/MgzDuXqV\nteyn7HmaUEr5DwH86DAM3z59fx8ADMPwU3va0VrWspY9lf1IE+7F4p7x09i+jbSWtaxlxWQ/bjpq\nhhqllDVruZa1HJAMw2C3jPcDDJ7BeB96yP0YowOR39qHrvdTdKhm0/v7sbhrBwBbUqbVVsjmdG5T\nysymNrU9LrNJfc/k81ai09Gp3tGprfdj3NEDxt3ILXqPeqEf67iZtM/lWVRHFS6/0yX6q5jPi2uD\n9Qwdot9an9mc7pe8H9vX127kW9Iz+wEGjwB4SynlDIBnAXwXgL+2D/1cJ8mGaDdD5xbUhnxm4446\napTAaMhOGAhcvxvyzgboDNS14eozMGi72vYWtvcV9RkYneh1qR5vQD8I8PdNae/1I3t+xcMwzEop\nfwvj7+mPAPilw7uT0DM8sej2ymOE11VR427pxmVdHe6n1lZmjHxejTba6xkTBg8GlKxNbTcD1jcg\nB9Ss7gYWAUH7vbFlX650GIZ/CuCf7kfb11d6PMS7pnc2gFjgNZBw3lOjgyi3JZ9jsTpD5bTAeV3t\nj9v4mqS801lFdXcGW4sCnE5ZGqTtcpmo/3Wih0ZWCkBZez2AlgHVXshDe9xeLq8f2NuxtDzEu5Lj\nIctGDVk/zotp2cwDAj7PVzD4WmlfAUdTAJciAN5gMyDJQKAVkaioru+W844/yVIaBhCXRoW4uVom\nKuqRNRisoPSSgWxELjyPdnRxZV5e+1HD0/JKFirhl5GUXF/77jXMzODdce3bAZWKmwM17NqS1mgp\nopJam9k43Himc+Nd0b6Lhpe1VEI9D7AdGHoiAZcWcBqR1XOk5FF4z1sLv0MHNTaOBLawPTLIAFSj\nLY0+WsvSkaPLSK3/+K4RBIN4D19z+GQNBqm0FlrmPVRqbLgLjTNCrLZrAGyPCrh/jgy0zxanoYbq\nZCbvvbsRy7D22fZqrX0WN7a1LVLQ5yyC6OExDo+swcDKMh4nWygqWTpRy0lB53jROsPQ71lq4ICA\nPa163NYefJZaZGWd6BjWQMLp2LPj4fRxKZMCZM9acBHg4ZM1GGyTnUyqLubMgHp4gtqCjfPsiWue\nkBd63GCkbTggaBm1E25vWfJP23BePIsKMtJR9cjSOgcIUebKErpHvehb+w0dudxqyRoMUqmF95lk\nN7FEezVP2TMV2bZbtKW5O9fJgOAK6t4dVNbtGigHomE316+J7u9HfceJOB5Gx1fHaQY/Bho5bcEb\nrRr6bkxnL+9L2TtZg8E22Q1jXAsXnefa7fBz/S0s/k+KcgWuPAPBThYnX6/b5Wgtes3VHSCEKOmn\ndyk6UtTdt8G6Z+f0eDZXoUPLQWQR4GoBwhoMrGRhfc2Al51YXkS93IESdaobh7WxcxCiQBAphtZT\nHVpLJACHx0avgQ05u9aMwa/t5TsgaOnb2npkfVU3ltBN77vQbeWWLqsDCK8DMKjladldgG5hZXf8\nRTs7ndieOmG4vPha23asl+rGKUWW67eMRseo976A2j0Oro/spi83Rztdzj27IGzgrFtGDLu1sszu\nyfWX1dVs15JNUquOI/2yPNJJr1fIPE307Rao7uVrfn0U241bSUZOC9x2ZIv8Y1DS1EDBMxMX0mdR\nTy30z8q4dh2foGPLfERtKzGTnrWyuoCwmlrtWmo5YI8X6AnZtY3WnWq1PDjTwRmIen7W4Si8h2U9\nog1NKbRvPa45PUcpbq/ehdDanxJ5LRDN7o/oAd/MCGs7KT0g53YnQqfMASnwu3avv9yAYLCTS3Je\nShdnbasrSzf0WGaoTgcVR/LxMfbSbqvTXYu7QSnTT49FH2z8LR24jx7AdrstmRFpBOTGP7uOGomq\nwKhcTyYaJXD5Gg9ycKBwA4JBS1okT43RVoPSBZSFz1kbNVHPp1tqtQjnGOmj+qosk9vWjMWBUk1U\n/5YRZHm4jo2mMTonrT7Z+GvEZM1RuC1Y1wZL6N2bau293IBgsExO1nMvQW3fPdvzdiSZCy+BfqPj\nc+zVWHgrMYsMtqQs65ClIzVhYK1t32m7y4BAtKVGxW3puLTIUE0PopwSyY7fYP2vUPnMEdTut9Ax\nV0L0+gHCDQgGwM5JmloYH7LMtlFGkumxZXgF3QXQ/jaweLch5/qXTR1g0ZgV4FoePvpQFt2x6i40\nb7Wr3lVTrowIhZRncUAQZTNA4LqsH49vtsOhPERvdOjq75/coGAA1AFBicRWOM9tZaRYiC4KNTSY\n76zHMhPvdM+igyvw3ooXew+pmglHCHod+j3L0XsIXt0BaIkjMzN+QHVxW6DAHICuYPEZkVsYUzRN\nmTKwX63dhhsYDID+3Le2ncQeqHeXQcswy8+GV8tDtf9MNHSONCVuOIp3R7Blv4RUqRmG0yVEQ/WW\n8Tog6IlSWtuNSghnxG8rNYg6W9gOBDOM4xnzuhvTakWL+yM3CBjoZewk9wW253y1PNix2npsMynH\n/fWkJiFKfLa2BFnXMH4m1rQNZ8xqID1j7SKWzOMvs52YnWNdszI67k4cEejmK0AgQCFEgdhJBkg9\nDmZ/5QYAg2ybZhn2PBPNV7n9Wm68E+kl1Fiv1jE+5whLl2KwOCCIdxfycv48w+I47ZU4MM7I2B4A\n4PMunVG+QH/YpXpxWY1Otsxn7V/buH5yyMEgC9uzcNttM2WEDnvyVpmjyfllvaB6ix52Wj9nnr42\n1Rouq8G5yIDbU35lma3JmmgEw327MXBEn4tcsrFQgpPXC/MCHNEpP9Qy6Bpfop/Xuwk7lB6SxoXy\nIZmhOqLP9ctRhGO845x6iKxvBQRu1+WkPUSgW1wKGD1bctmycduG3OaykVprq0+NP2urx6B4Xlxk\no3PpHJGOH88d97GMgV8/QNiXv2Rfy1rWcvjkEEcGyzD7XKfFFteEvYMjiSJMdGy4eoUs7eDwM9v2\nBLZ7R/3eEx04pj4bG97Hz/RqbQ1qhOTqs7g0YSe7Gq59jlb0PoXauLM+GeHqdMnmP+NpuP71iQ4O\nMRj0Sm2LRsm/3m08rV9jj/UGF93ey/RVycLlVsju2s4WVzDhuhiDB1BuINNXDcWlD04vdz0t0Hbg\nU8vHM3EA5VKmGmGpdTJ9nLPQtpST2X9AOKRg4KICN1C8KDMjygY52wJU76ELw3kG3o5SyTyBE7dg\nHHmasd2a73P0keXiW/RSvsORXio1UAgdgO3PaKyJm5tMr5p+PeAf/EwvYVkbB9YziwwUEFuR3t7J\nIQSDWnqwLEG1ExTPFlBEBxpGur1q13fvFmWW5mTbfVrXtZd9V/KVAY2NrmXEtfSIWXgHjK4tFrfV\np2RrSM9aqAFERrDWdpxcew4MuT0+vv8g4Ho/BOLU7Z0AldrCaHkz1094EE4blP1vgZUadJaXO8m2\nsLSdMOIsfeJ7BbaQG5vrl9vVa1GvqFuYGrH1GJimX3q9Pfc7uEhJ9YtyDFg9obtuVTIvwXrW6rbW\nci0NXk4OERi0UDMbjMwTs9Ta07CullawZ9cQks9zaK5lVbIcWsNH3QtnqS1AlxKEEXGKw+e5TY0O\nsmiL9XUhN5d1ejvQckCgZdyWr+s75pWNnfvWqK2Xy2hJFulmaarW4++7A4RDAgY7AQJnEJkh1Lyq\ny72X8djh9SJaYE5A88PaZLYWHxtsJi3ykndLnGeqjb1b1OyZa8x7FhLXUrYess/xKVq+FSVByjkA\nz3aQMnHAye3xLc7LkKG7A4RDAAat1GCZBRPSE56pMXComeWjG+YYg9BR9A05e/uMqOu5bl6k6iVV\n/56tz2jT9euiHgcEblekdT0sLnJwzkLHOSOY492lNdwHR3chnBr2cA4KlDUCmPt29eP73nEKKw4G\ntYigtWiz0NF5oZgA3kKrAUlG4Km4xeok8youzO0No0MUEEIfFzm5iKC2Y6C61rxsr7Ry5dqcZtxA\nK8rTnRSXOrHo9muPR2bQUCDka4i+9PH1GRG6d4CwwmDQwwIv4z14sHXC47i+tP+YQCUI1ZMvKzVv\nG9JDHGb8SC0f5WOZIfZ4JPWMzGHsFhhYR8BHZPx5mVTSRXOcckU9fjCMSxcUBLM8PwMuByg9fNfe\nyQqDQU2yhauIDvPdea6eqEBTBA4te4zNHe/ZissMs0YSZu1rW+p9ssWb6alpkOrDDHovK65tsG56\njQ4UsjRBCVfmbrJoS79fnurEg2J6jdRFkVlkEPqovj1ywxKImitnZbLvbqG6vM4tWo0KGP1boWjU\nqZGaLnfUvFNJzBm2RyQ1yUi0aL8Wjmee1nkvt6CV5deULjOOOMe6aYjfYxwZICjJqlELRzIZIPP5\nLG3KjDIjWVt61+R1wxn0iAvn+Xst1HIL1skmvCdwJJPLx7U86xNT4J6SU/OETscsdM8ippqRRl6c\ngUotT64BAffN4+pAR9M1R4a6tllHfo9+WfcrmHv7OF8zsJijGdr/0lwbIxe5xJjUItPQoTca7ZcV\nBwM2rp3mTfq7gYz15YWnuVwWKjtPkKUtrkwsgNZvG0Inl5JwqOvGqEYAql5RxuW1PUsliwiANtG1\nicW7OCMUd6mMAm4NgOMash0DBRx3zIEKE30OmOJ6eVelJa4s61xLg3cvKwoGLePPJp/rOY/mcmhd\ntFrGsb8KAplHbon2W0P3uLaeHQp3jhdnq54DAV3cKlm64aKhWt8M1qzHlaQc981ePyP4gMXnQsY7\nA0AWwdWiw5r0kOFZNNCas72TFQSDFuI5LiALp2skTQyyphkaHbQmrxehlaxSD8XnORpgUdKr5vWz\nHQp3PosoHGmaMfxZSsBRj6ZTcSweGMvt9wCsI1OZX4mxZIOOqCM+8+6R44uinktnYh1qpMf1s2iT\n9VVZdm3tjawYGGQXv5MUQScO8ISP5mhqAI4wbPEL3FcNTJTAYr10wWh6E3rXQCCbXvaWSqBl3IoC\ngaYDPXyF62+D3hUQs50USNnQTUP7eDHxynMeL/79hYbkzJ24dNHNTUtnltpOhmtz/0x2hcCghYLq\nWZcVZ9QaciuKu7TC1ec6XC9rh/PNVr81D9/LRGchvIa9WV2OoBQAwpC0bWDxjsuMmGOjvIicPKyN\nSZa21Dwsg48DkygTXIZbP73iHIOCaE1XTYncPLXSzLbsCgxKKWcBXABwFcDWMAwPl1JOAfh1AF8G\n4CyA7xyG4fz1US8j+pZtg/tuhdu1OjVdop3LWFwIGtK6OpmR6IJtee2ajmEcGjlpRFALrV1/3I77\nrP0so7NrK4u+uK0YN+0vgMCBgduVydqP9pxDcztSXI/P7RSM+mS3z0AcAHzzMAzvGobh4enY+wB8\nbBiGtwL4F9P3DmEv5S6aF55bgDzYNfaf80Ng+6ID6kCQhbCaV/aAklu4mRdrCYexrv9MZ2D7+KmO\nWV6dGTFLawFHm/xfBG5OQl8eZw3f9bj2EW1rNKNtHcUiCBxNXq3Ij9t385HxL/FZ15ojOPdO9iJN\nKPL9OwB80/T5AwB+GykgOKR0IZHb+uF6PYShG1AFhRqplnlZHcJlPIjTj+s4PoD1VABVcovbztIL\np6eK042JOVdfCcFWu9GG6tYiRjPyU+u6NEB3c46RDgwC8dmNrX5noGzNr85LZo4BCGELuk4zfVpj\nuczZtgwAPl5KuQrg54dh+AUAp4dhODedPwfgdLtr/qwLxOVSOok1ss/VUxBQUVDhshkYOMKSRftz\nizoLndnTa3/LhOMu144FlhFZ2TjVmHMF6dYyc3nwTpamA/MMjGqEp4JBAAEDr/PSjoCcJXUyvdz6\njbF249TiGELa47lbMPiGYRi+UEp5I4CPlVIe55PDMAyllGG5JnmRtcIiXmwasjnPr+2qoWxhcSss\nY4ldiMnlQ2qgE4Y8w+LfnjnAifKOm6iRna7vHm+pht4T+mckHksLUFTHZXLkaKsVWTq9QhQQ9OXy\n+GhPAZVBIAv1M65FgazmCFvSb+K7AoNhGL4wvf9JKeU3ADwM4Fwp5e5hGJ4rpbwJwPO+9i/T53cD\n+NrpMw9AD8Oq+V5W3jHULEreRZsZm66EFgNTxm2o7htSxuWMvC3GesF85mvpzS+djlko2hOi18ZY\n64bXc+3VgCBr34HLsvl1xptk+jAPwZKtxdYOh7sbNbtejRT0+jcAfGp6tWXHYFBKOQbgyDAMF0sp\ntwL4NgA/DuAjAL4HwN+d3n/Tt/B99Lln4jVXVsIs81Zcv2eRhrB3iLYcmtfyPPUWSnZxG3EjzGUs\n9uPKZiAA6VNTEm5LI6Ysnciinp66LUPs8fyOV2qVb0U1WaqmEmDlogsmIpWcZZ1V/1Z64vTS9mdJ\nOe475OHpFfIrpo6vuYycBvAbpZRo51eHYfhoKeURAB8upXw/pq3FvuYcsecWmRqE/qCmxTm4vEv7\nCIly/CMi5x0VNPTvuh2v4UJONdCMnKxxEyEt4+Ex5rHWiIffa23X+qsZQFY/wJ9zcGcQGllwju3A\nv+blnc76A6bo4wq981OjnSFrmzsRviaXCoXsJKKa97AjGYbhcwAeMsdfAvCty7XWE9I6Y3K5HLep\n0QAPXG3QHA/At7eCjimJx/1y6Ki5p4aRM3lXyfgRVz4LQfm8AkGLF3AEZ+2Go5BWfs71nAErODmw\nV69ay7Ed2HFZp29cZ0RwcSz7L4xMVDd01ufIRCO73aZWi72sgNTIlcyjtkJ0Desz0s+JEokhyg1o\nVBD9uIXC+SADAeuZEX16TdFGi6DLQntus5U+qQFFHX2EehaxOAN33hzmPc5HOy4/dsfc2Lj+1TiV\nvOaUS3kkTQFr/YT+0We0z1IbI66XRWW1KK4PsA4QDNxC1mNZNBDHlpUakKiHiXCTF3sYI+uQsfN8\nPS6v3zR1WLL8sZXuaFjtFh+DQM9C0YhAASzztky+OgNXYlb1574doenm00UPNfLTjQ/rxPq4nyvr\nHLs0kK+Hy+8kSmCpcQe18l5WJDKoGUOWJ/NCi3JOskXgzoUusaB48rmf1iRkg84LmL2LLjpe5Lyw\n2Fs5UkslC6szr+bSruiLjTUDEhdxuTSNxzhLdZz0lmuJOpQAzYxj4PdaWtTjqLLrXSbK6Cm/bPsr\nAQbZgs68uCJ4xqiq1w5xjLy27YiqOJ8RQuode4BCw3QHVC7UdORZVkf71mvQPt3YKL+gYBjlM+/s\nDK2V4vRIdr09HJRKr3685vR9L8ypBjbsFDSlcXWWA4wVAANeKC2UXWaSs7zRRRrssdnrqdG1GPNa\nPs/e312H5nw9ZBATlFkOzbplQKB9ZmNeq69RG5eZwV9PbX4yPVqiY5x5+01Ttje/VkDQ8csITBVd\n+wr2oHP6eQu5+c6kXJ+sABgAdYbXse4uh1WPFoSdM2RdvBkD3RuOO33d0LZy9Mw4oy7r5n7+zIsg\nvisQZISfA+JszDXqcnqy9Cwztway0D1ru5fYc5KRzb1994KdinNI4TiyNLhm6Ap0XOdQcAaZN2xd\nQGasGbMa4nJGlRq76/p1Ib4uzuxaliFIs7RFiUPVyXk959VcfwoEWQrjIp1a25keGZBm7bfadpJ5\n714DzqIibceNedTr0avGiwVY9hKI9TFbATBwi8uFT5mX0AFTD5ldokPaDBR628wky+8cz5AZJzPo\nNW+Z7SBoW7WUrBXaqgE4HidLATT1ivb0J8EtAAWdz5xGTxidSW3s3Bi4estElDo2vdFpKxJdToMD\nErclVQvd1YiynMkZ737ITiYhC6EdGNYWsnph1olD+IyPcXqpAfeClquTEaKqyya2/0SY23Vz3MPE\nxxhkurXAvebVlWTNUjk31600MAv/M0CstefkUKUJTrLcp0WwtRY3LwZGZLeQWcea99ac0y0mbV8X\nlYsIuK0NbDd41nET242QxY1J7bq5bdU3i9S4vPbFemXPC+A6nHL1iLt2Baks7XH6uutRQHA7RDW9\nnLi1wf3rZ21XHYgjsOuyImDgvLgSXsuEgVy/V2oTkQGBG/wWiZVFFD05pGuHv6vXZr1qHraWd2ZG\n0KNntNuaw6MYHyyi4OwiINaH9c44CtWLCecoyxFC3DmagV0Asd6hmu3Y6Hip/q59Ph9gU4u69Npr\nUUYuKwIGgEcyHtgsbOXvSr7UBsLlYUy81bZusjDesfasnyOAalGB0y2Eo4NaH600xm0D1nTl+tni\nd31l+nPoXROerywF3E3+rGN/MSnDkZdLZVskZy8J6lIdfl5mjP0xbI+UMyCrywGDQU3hDGmBPN/j\nQagtwp6Bct6ltsPA0kL8FqMMbNc7izayvDG7xloomUUIHKr3kIwtcX2ogWdA7Yi1zEs64fo10N/A\n/DFofCzTcSfjkUU92fgEWGfpTi06ykjW7b0fsChjytFALURqhV8Z2eT64jaUyHIgo+IITtWx9V0X\ndouDaOmh19+TtwJ9IKW6ZfwJiy7WLKTlec/SD9dmj7BRqJHourpC5zXcz/ioml5ubrJIzh139WOc\n9BeykO+1SNf3eADiyL44nn1XMi0kIxtVasblJiHzSurN9LPjHFQyAFA9a7mzCter5ZY9099DrvUI\nG5Maf/zKU0UNjiPFMI7WOGQRVSvCaRHUqiOvSW2jVj7TK/ueOYiwIQW2dkSQ9XQAUstvap49Ptcm\nM1ssbpHF52PT52xoLtPnIL0492eGeSdbPbpIaxyEayszVPZ68fko/EKtSS39qrWh/AgDQi3iydpj\nQK557B5eIwPP2m7RTnYMdI5c2snX43aCjmFcg45PUOeoutUjywMEg4zw6DHuHja7dWn8JCKNOmKg\nlTHWCOEyxslx3rfmFfi4i4r4fI0IdYvVeR4m63T8jlK9TWw3HCe9EYELhV0k12LgW+JAPwOCrF1X\nJvTgseR1w5KF91nE61LfHh5jA9vXDEdYO5cDBAO3KFoIrRPJ3igbvIwIzJAyzm3KZw2/dDIYlZVv\n4DJ8vMcr6vna9h7r48QBAoNfi5Ng43akbq1/p7Mz2IwEWzbKylKElrE6DkTnykVqnHapI+jZ7tuS\nchkXozs+Ou6bWHwyU/8YrkBkkHEDSsCpqi6E4nJZPeYoal5ZQ64N+Meo8+deDoH1XSbvZkNuLepW\nO/zeQ0wC9cUMOtfj5Rn4nUdzS9O149pXIHDOhUPyWjTUInF1LjKHVWubv2uqsCygRbSqOrblgCMD\nveha2RrTnU1k7fJqhhptMjpznp21WyMkOdrQ8m6B8WdezD2RgfbPYMXt1dKPGrEZ7boQuMUZhB4B\nAm6hZ6QntxX9ZXwB6xTv+o/KKrXopjfdUD0dsZj1UTvX4oxcWR2nleUMons2OEd6qSd0CzS7DCUb\nud8e4dAuIoWsTefpoxwDiO5OtJj/MGCu7yIDlloKxDqocS4jWS7Mnx2Az5A/TLQG+DDnMtLScUGu\nHx0nF3aHZOvHEaLah2ubgTmTGrmbkZs7lxXYTQip5cYZkrfRLp9cBhdHFGXi/ma8xlhn6YRKBmyu\nHW2Py2q6VEtH1GO02O6WgTkQ0Hw6IoILVOaYlOc2svSPr4GvW4FA82sXaWlbXI/PsbSAK/qopY0q\nrWjFSZYyTQUfAAAgAElEQVQ+11KXvpauo2RMsi4k5wmz/I3FLWQW9vS1wXJbPTXgWDaccyRiBn61\n6YpHuTtPpOPlyKWauPqqT41cDAlyi69ZI54QBRYWR8T2Sk/ZFhmqa7TVdiuV0TKuXEtvjZ6Xmd9R\nDjgy6GG9+fHiLhStSTapkZIoSah9KAhpKM/He0PtDACiPxZnbCpcJsZKvSWw+GQk5WmyFEb15j57\nc2bnuTZJV75XI2QnHJDr20UNmdSIQgd+jpjeiWhqmK2vmjgQ4LZqXM72Vg5AlO3lAY/F0sMUg+px\nmWhXjzlhUMjCrgCPTG+toyCiujNfwP1wWK3g5PrS8JmJStbZSYt7YN25f52bALeetjYAnMIicNVI\n5BZhrPMRdThq0LnIUjo3vi0j2k1U0uJGmE+rrQN1nm5edb35Hg9IsrwuLkTTg1qe7chH7qfnGC9u\nN5DA4l9tcf9RR6MMXYQMBvqXXHzdmjtnKY96fxdZcB+1qCeLbHSeNuTFCzauU+eMjZQjg+AP3Hi0\nog8FWJXMc3Pk01r+NePpqa9919p1a8+lXNn1bmJxbHkO2lGB9nxAUhsYPp8RdXGuJ3dz5xxhlxE5\nDsC0jiNutF7NU2vYmF2rTrKy2VumDOvIOl/A/Ce74V00ZYjbsOMns9kiU1BwLLhGYS9N7cf/VCjp\nyOLGsUUyKrfQk6Jxe06HzAE5R8KGqeCs88268BjW0mSdV57bPiBwmq+AZOF1fHeLW3Pw2gBki8Ll\nXRqCOm5Aw3qV2uLlNnQiXdrCurjQUReeu3VWvXN45mD3T1D/8au94/T9MrXT8t49CzF0YXBhrihr\n1zkFTouU9NV1o9FKfAZ9rhGJGVHMc5CBBhs4H8siU9Vfr30L9fFSfb28oVliLWtZy+tCViAycGE1\n4NMDJGWdF68x873iiK3Ic6NN9kAaUWS6s2fgiMRtWdUQ3Xkcl07FeCgn467pBJWbYfQ4d0znngZw\nFmMUEc8vBBYfFRbtOR5IhXd1gkNwxGK046494wtcdMB1VSenJ/cRc+YiLdXVRQdORxY39youOon1\nlqXPHDXXo4MDBAPXtRpBDL4j2zJykQdAc/AWKZcBSITcoc+WnA/j4lCttlhVL25Lj9dSDz2mAJDx\nBc44XFqyAeA0gDsoan4Gc35BwYDD4kgxeFx410H1OYE+o2Fw1TSJCdPM2OKzA0gXpke72Xy61ISP\nqzgiN0Tnxs1THO/Z5ZjJqy4HCAZu4vXYFkYv5PJ0Vz7LJYHti92RW9FnbViyQc1IR36veUknmY4u\nggJykjX6OyrHw2i5zVg4YaAnqJsw7KNUNtOLz/MOiivvSEsnNfLV6RN9aPTJc9wChDBq3UnqJ+a2\nA0VGHmZg5HRyoKSk8Za86rJikUEILx4mwhgRNdTNyB434MrQKurHudaEZ+d5AhwQ9KQAtfZrkqUm\noYNGU7HIAyjiVuGjuAYE2MLI+L+IRYLxJSwSkhmpdxFz8nFT6gAeSLl+dsyx7q6Oi3r0nB5zukW0\nU4sOuE4WHXAZPe++MyCEDrVUgueYwaG+ng6YM8g8AYfjLVRzqNnqM/NgIRrK9khtu2pZL8JhZ+bt\nge3Tl7HXR+mlue8Wxt2C6Osi5ov4+FT/mem1gfGGIWAEggvStnt8WYzLZTkX16i/XKyNeZYCxLna\nelp2PlkUAFrz6niCGuelbdbWcgs4A9CZS8ii5UVZEQLRkYWxSJwBsGfRiXLI7STjDVQnLs+e1SFw\nNlG6EDkCcWUVyd0v/NyfjrD+0c8mRqM+g9GQrwA4B+B5On8KIzcww0gSXpyOF+rrOOZkIqdvcT4i\nBo2yIgXYnPrmslcwgkoQlQE+fC0sHAk4zqUXdF1k1pN+ZNEir1Md+wzAM2KR+3PphZbRKHBLjvc7\nohUAA5dXakQQgxCLSg2TWeme8J77c4ugZtg1Ik6N3A1vRkRl3obHJMopyRfiCMyjGI39DuBOAOeP\nAbMTGI1wazp/GvNfDt6HESyOzpuf3YH5LkPBPNyPiOAYFg2eJfoPnTYwTz8iumBvqLyGjkWMAV/7\nMkCQtRkS6yJL51w05gDckYp63omuoSyqYUfEemZpb1sOEAw0x49IIP4oIsu3aw8XCSDQnQdgu7Fp\nHzpxDm0VaNgwnXErGanAF5KBhgJk9H+U6rFhaGQj7V6KImHIUmYDwIy5ggGYRXSgHuo05nPBqQen\nF9lizIC4FcbXdgiuXYD0vYzomnR6Oz0yh+Ki1Fquz3V79XdAtDNgPEAw0Jw4gIC38IDFxR4PH2Uj\nj3P8OibnVHoXYouEy8q7FITLZACk+mmkxHpmKQbzAqHnRQCPAS/z35iFwc8weuhNYMZ6x/Hw/gEK\nocNpAHdhzhNcoHPxQ6Q4d47qanpzQt7jOpYRBiQdewdIrbl2ZLWTFvfUKu8kC+9raUWc313fB5wm\nuJydB8Mxq5kxunyq1mdmTLXtMe2H67k+MuGF63YZtB1+1+vMwkfuh397cArAvZhHFxemV7TNIWeA\n8wksRmSF3uO3Chy9HJvqzDDuQLw09R8RIEd5fJNTXJPb9dBrU+EdjRAGB04leRxbIOKkx2w4he2t\nz2vfgVOPZE6u3cYBpwm1QYl3Fwk40XA8A5qsbrYr0cq7aryAE01bGAycro7FjnBcgXFD3hnAOLUI\nriDu+ruIa9EBTmFuoNHGFtUPIBjoe4DFOcw5g1OYRxBbGInJZzAHBQVBdQD67sStiRaYc5u90aCL\nypSo43ZZl1ZU4dZKjfir7TRo267dXA44TXCeP957wAJJmSyPd+I8M5NCLcIv08317QDARQMKQLW+\nMraa5TjGkD7Sp1hssTuwifk9A8zqR8p2HEABbpkOvxznLo7t3lLGFGN2AmOUMYHUSYyk5WwTePrN\nwOwYgCcw31Xga87+HszNdW+4rEZbWytZG+54i5fI5rYlCgIKQtqfAw11Nq6elxXYTQDmHo492Kac\nC9HQWCeZw71a7h5tZw+DYKN0hF8t7HLRBANdXKcCkEN75UbimJZ1v1pjUu8YRqPmPHgycpzA6MkD\nIML7T+BxC4DbMBr3BoDzAJ4L/oa5hACXjfE4By8nAbxwF0awiJ0MHv/QV/XX77UlmxGTjljjOVfy\n1bXLbWW8j4K8AkItyuQ2XZqQ6Rbi1om2X5cmGJRSfhnAXwDw/DAM75yOnQLw6wC+DOMvV75zGIbz\n07kfAfB9AK4C+KFhGD5a70FzOg5JQ0VNE2oXlp3jCdd29WEQ3I4uMM49a5IxwmwhsWAC4DQScWGn\n9qFtRzleeHEjUTxT/ziuGfptmDz+sdHjhwobmAMAv4e6LwB4YXN8X9ilOD22e8t0/BK1iYIRdIKn\nCM4g0owACb42x9UoEGvIrhFBr2TOgwGb21cuh8ffpX8sWXRR69+JA7L9SxN+BcDPAfggHXsfgI8N\nw/DTpZQfnr6/r5TyIIDvAvAgRpbq46WUtw7D8Nr2ZjNiRQ1fBzlE0TrLpTT013bZMLntGb1i4hm4\nnGSIr2VceK/pQauvTDTMjBD8wvR+F67t+8dw3ILR2IE5IADzMP++6fMt03lMx84DeBzAkxiNfoPK\nzabzLw/Y/sMmHgOOkCLiYv35s5KAfM09INAyUOeBmc/gclmI3uIxXBoQxx0/FFIzVbeuXX91aYLB\nMAy/U0o5I4e/A8A3TZ8/AOC3MQLCewF8aBiGLQBnSylPAngYwO/6rlsDl5XhcMrxA62QyZXjthkA\nMg9QQ98aKCjxxx7CLQjVv3dy4zpie4/BZWNeJDz3BkbDvxtzbmBj+v7A9D4D8BxGI79Fyr2AESg2\npuPx/nTBeO/CBXrF/fWg75E+8LU7hl9TJBcVqGQeM2tHPTx7X/7BkpKSTOxyf5pKcN81gtMJt9Wz\nU9HPBOyUMzg9DENsHp/DtdgQ92DR8J/GGCEYyYxJB1LzamBxADNiTw1J0bKWV3FkwCw372zsxkj5\nv/Cy+q2pqZ3XMDZSoQjHQ6YbvF4+Cjy3OXr1uzF6/Vj/dwN4B+aRQKz1l7EYOTyHMUJ4AYsRwkkA\nZzeB80cxZpR8+/EGfQ9Q4LtMnag3jWNZWqGGnbXnvLyCfnzX8eX6m9j+Gwz9t+tMFIgyXV3q5NpY\nbk3tmkAchmEopQy1IvWua0akk8GeorYA+HNPmM0o7Z5bwG26fF8niXkA0DlON65g+0LjPlh3XbC1\nBc5eV/fx9QaiYxinpyxyaScx5wjuBPDAFm678zwunzyG13DreI6n7fxUNlKIGeakI6Y2zvOvIvm3\nDjrmasQZGagRAYtzJFmqWQMRF3Lz/Lu1G7dqR/m4+SpzIA6MFIiUXAz9lKvQMiE9QLRzMDhXSrl7\nGIbnSilvwvirF2DcSL6fyt03HTPyixg5xtcAfDWAd9I5HQwl/wBvhCwcrgGLA6YAoamAtscLIzxX\n9nft2p/2oyQhsBjCZyDmPJduzzGIxW8O+G5AuvYw1lvK3IMzb3AngHcAt595DidvPo+b8Apw69jk\nVWzgKo7gMo7hT565C3jylvllncQIDpewyD+gYOQrQvgHS8cxLiFngAqIeq0qmTG7OXefM4n1FADm\nQMmtUWD7WlF+KI6x7rHOnJ7KbWiUyXV+f3oBwBHUZKdg8BEA3wPg707vv0nHf62U8rMY04O3APiE\nb+K/mt6D4GLSDFgMu2p5dRbmuclxEgPpymqkEa/WwyddfWCeUzoEz4jOWh6sYMd5bYznsYn4kygg\njP9O+hx5/i0A7gPe8rbfx5txFkdwFQBwEudxGudwE17FZRzFF3APPnfvGfzRnQ/giyfvHtt6GvOn\no52HTEPclchPQc4eatKauyxais89QNBPrs2F5z/60fbc2suegtUCIgVD5wDVRqKPDYyO9qupnQ+k\nPTXBoJTyIYxk4Z2llKcA/M8AfgrAh0sp349paxEAhmF4rJTyYQCPTdr8wDAMtRSCLiYuQJlXNdbW\n4DHrr5fX2gVgw3PhZfYHqjUP5dIcruumYGY+K3kZCyAWWfABFzAGY1cw3gp8FzCbACB+dBRsP+8A\n3Dm9bgPwAPDGr/88vh6fwAN4EkdwFTfhFdyDL+AMPgcAOIs344/wFXgTnsWDN38WZ7/6DB594CFc\n+t03zi/phek9uIZZXEeMQRCKnNbEzUzxPVuibo55fDJPWZOM9WcJrx2fs4iVQUKBWyPETI/MSUUb\nNXHRVX8NK8Mw/LXk1Lcm5X8SwE929W7VcQaWXXjNmwLtwc50cJ81dcnY4IycZHDgkJDPOUIxS4H0\nc9w7cBxzI3t+bHt2L7bdV/Ay5qF8gMAZAA8Atz/0HL4ev4d34xE8iMdwEudxEufxwEtPo/wxgNuB\n0295Hifxp7iKDRzBDE/hfpy89Tz++UN/Hls4MfYRZOIljNHCpbiuiA7iCUhsYBEpxnXH3ZAaMrsx\nrhGyNXEpWcwjczvsCHoe06Zt8ueaMavDYD1deqzrSSOiHn5uDwjE3YsSOvxeQ+caMxyfey4v255p\nISunLE5PnaBoJ2O1tW4mrTbijsKXplc8QOReAPfN13hEB5EWnBxPbz5wAV9x85N4AE/ibfhDvBOf\nwd1PfRH4PIDPAngKwCngLX/madzz0LN45eabcBUbuBMv4lncg2fveBP+7Ve+G9f2HV/AfMvyyTum\nTl/E9udWOELWEYA6Fm7sWuNXE5d+MMDre3zO9MvC/BogzLAISLy+1OkpILDTURKyvu4OEAxaYU6N\nJdbPtTZqEoxstKdttthf5TAci+sWs2uX6/BxRvlapBJeFhi91r3T92cwuuWJxX/5jkXO4CSAr8S1\nG4tO3jFGAXfhebwJz+Luz30ReAQjEDyOERRuBvAEcOujr+HWr3wZeDtw/q4X8QD+CO/Co7h473E8\n8dxXj93egmlHYrqcs/EchBcnnbMHsmxi+1YoXzuLy9mdMAMf37VvYPuj2DhyATx57PgfFd5ybEUG\nql+NAFWw0pS2p88DBQO+ADWGDAgcIaR8Q3bBPMButyBbJKyvRgOunSjr2m2FbXoNCgTxvkHfo03+\nezJ+NkA8SYgW+CXMbzE+M73uBHDbgGO4jGO4giOY4Wa8OtrsExgB4XFga9o32nwKYxbyMoDTwE13\nvYI78ALux1O4Ey/iCXaod9LlzQA8HT9o4mcvRhQD+HsiQnQNZDtKmSOpzTHn8ux1nZNg3iaLDFha\nnBevY7dG1Vlxey4t2EAOSl5WJDKobc8oCGRetZVDavn47MI+V553PVxq0OoL2K5/WIwSi1n+5yY4\nxi5+OBQPE4nXJsYtxi/Dta292Pq7BWP4jvHz5smLOInzOI6L2MBVXMaxscuXADwOPPHE+PE4gDMz\n4NgpAG8G8Hng9P1/grO3XsIRXMVVHJkP7bXfPmC8MWkDGHc12NC2MEYJ8bi1ezG/UUc9tJJrCgo6\nP73rguvGS9uqeeea4WUcgNaLMgoILkJQoFRyuYeEXJQDjgw4J3Leli/Q7dW6i3TRg0orp3SRCv9v\nodOx1Yfz+M5bZPprH3GbMe8mnMJormcx3g4Yv0V4O4AH6bcFL45dPz49F/HMePzmW165RhYewVW8\nipuu/eDowufHTOECxozivlD9JQCfA275UuDmr3wFl6fY4hoQMFH5MkYguvawE57TID31GvkPYfkx\n7SyZ180iqzjm6nCqwVFHGOdMyqhkJpVFKlq+ZpJ6jlOBbCeln0hdgcjAESoaWvOPWRQEHDq3kD8+\ntyTKZ3clunSllabUjqn+zhNlBBkv2JDjmN94RM3yQ0fDQC8Bs9kR3Iw5IBzHReCLAC4AL74ysg/x\nCJQrrwAnrgL49xhTiRdxLSo4gqvAbcN4Q9PLGCOC5zD/heM1YOUHrGRgXyPZgHw+W0aSifJIG5g/\nKs4BevTd4sAycM9IUF4LmWPj9CRLo6O99no/4MiAxeV+nB9nN22o9ORKWZ6pZRgMNPRkVr+3v5Ba\niuFCYRb+4RE/SfglOn8K8zs6Y+/+AnBp4hJuOTZ67bsxvs8AXAJePn8cr9x+M45gNub/V58CnsVk\nwPMrib9TOc0/cpou4aYJTG67+wVcOvvGEQT+AHMwAEh3/fky/whIQ/PM+FpEczY/2doBFm/pvgPj\neG5gPvYtfsK1q+mfA4faNQOLa6H3GltcxWIrByQtIkUvzpFt/N4Kt2IyWqyvLjL3t+ZRVz1OC32z\nnJG5i5CMg2AykD1D6Bi3IQdQBLF4FuOivm/+46K4v2ADUzi/iYu4DVdwDBu4iqOXtsbjtwNvvh94\n+KmxlRMATh8Zj+PU9LoLuIjjuIoNHMNl3HbrRVzaeOMIAE9iBIO4rFkY16kpdQleI36sxKQY3+AT\n1xljkfEEXN6Rbtn5EH7AznHMH+EWW6L9ofd2qUUo7vgykQOvw9ZOV1/v10lcuO0uXElF9QS9l8D0\ndiaOw9B+OCLgRcMLuCYtEsmVyzwfj0388IdD2vhB0Ev0OjFGCPy7gbjXYAacx5eM9wvgHjx/+ync\n+5aXgG8EcAfwtZ8FvvbZqc6bAbx1ev9y4MKXbo51cBcu8p+hUBoyx67NscHY0bgE4LwDVBcZsadz\nn6NeNnYsSh7zuLK4NnitaJ+urPbZG61kZRgQsr6WSYkPFAxUeoi4lroZgeS+ZyxrzViZydXQjPPJ\n2iS1GGgmVl1dBiHOte/A/JfkFzD/X8Rw3RHgnwPOHwfOlu3r+WXg2S++CX94+1txB17ASZzHO9/9\nGTzw5U+jPIsxZYifpE3RAG4HXr4H+PSRd+IxPIjH8CCewv148cU75+mFDns87yB+5nw+rlGjMHUC\nTnryYUf6hjKc8vF4Mji9RJ+1XRVeB7UytRSgx2nVxI3hSqcJIVnu7S7IGWCU5XocMrttGc7hMtSv\nET5ZXu+2iVy7Wc7Ji6CWRmV6c3u8Zx/14vNF4PwJ8soYbxAC8PLsFB57x4M4dvsV3IxXcRHHce7U\nk7jz1Is4/o6LuAmvAgBexU24jGN4FTfhHE7jM3gnHsVD+MOrb8NLZ+8BXijzB57EcEQEEsei70vA\n/Ke/YbgXsX3ueHx6Q3UNmV0+zhFC7NBwhMgcgks7Mm9f44Zq+mrd2s5HjImmwCqtNbMSYADkbG/G\n0rbCIxZFevbybv9Xj9VAwPWbRQusZy1ErO2OOCCMc/FMwdA9OAP+A5Mg607M9/4vYf5QkjsBfOUI\nCI+++yHgVuAyjuEijl/7teJVHMEruAkA8CpuxmUcxbO4B5/Fg/gUHsJLj947Ast5zLmCSBEiIoh+\ncRnzR6fPMH8+YtxNyQ9LzUL4GAMXrisQ1Dw2j63+/Xpvjt+KCLQfLpdtW9cMXPvLHF+crxPdBwwG\n4eWARUPT7Zte1j++q4fWPvWzGm1P7pZNki5GTR80+tH2uCxPcItc4icFhRHF7xKAa/+3GHcnhv28\ngPnNQPHQkluAS3gjfu8dX4/zd5zEi7gDZ3B23GoE8ApuwpUJJM7jJJ7FPXgSX4HnPvnl487B01Ob\nT4N2ELAYEeDCVCD00z9qyf7NWQnbZfiaXslI3V7piVp0vjOyr0eXLJVuRSqLsgKRQbYtAmxHNz4G\n5KifIbQLDfm7y69iIpbdR9Y23VBzKBr9tjgFrc968lhsYv7nKGFkROxde6IxFrmDs7h2o9HWcyfw\nB/d9HR574EGcOX0WZ/A53IMv4CqO4Fncg6dwP85fPYnzL5zEa8/dOv52IR6QehbTA1Hn7Y0gEGH3\nFVEgbjiKf4wOjiNuOtL/Wojrb3u8XGb0QqOd2vxnhGNr+y9LUfVzlrK6dloOJJcViAxYHMI5FTmX\ncwRSxibz4PKgZXm8GlcmG0mZFnmk4KOA1lrkmccIMi7+Ni0iBf43pCuTp6bHoEWk8Pj0fhbAbcBr\nt92KPz7zVfjjb34r3v5ln8ERXMUfvvg2bD1yYv6sgkuY3/j4B9M7ICThOYzEZvAXscV4EXMwiL9/\nC7mIRamBuDufSUQT17Y4GmWdOOZeJeMTNApw9wRomlDbiciktlO3KCsQGdRE83QGAHcDiCMiM6DQ\nQWqRkWGcLjerhXEt/UJHl1L0SMsQ4hquYE7SBUEX13QcmB2bs/9xC3H8fgGYnk2wic/O3gVszIBH\nN4FHMff8M8zTg3jCESCPP4t7HkLXAIUAL35315SlWVkqGOPNoY86B+YHdsLgt5wR68N1os1eTkzP\nBYjViHDWp532rBAY6ALYlPcQXjQMBsyaZ4OZ7bvGYHLYquSTCvMZGdvLbdZyfued3PXXFg4bRrzi\n7sP48dIGFn/RKNHVpePApempSHdjziOEZ38UwHPTX6nFLcZsZ3yX4W2YA0E8hRknMN4HEaDAvAD/\npdtnSb+4vmUfNVdLIVxaxqy8K+skI/3ivZYCa3lHZodkZGd2bW7ttE39gMGgJz/WkJ9Dd50oBRHt\nQ/vVMi5UU91UT2a5e7iFTBzp6XZYVDSS4TGIZw5GmQi5MxJzystnJxY5hUgDHsc8GlBVX8YIBOen\n77GFOIt/YY60IG49Zp5gE3OQmu6FwBZGwjNAgv+AhaUWGWXGpfX4jk72oL1MvoozRkcYap0tU74l\nO00lt8sKgQFQ3zJkL+bqbWL0Ovqjpt7+MzKSDbMX7VX3Wr8aOWzK51Z916c+jpzb59CY+4u+pj31\nl4/N04Tz0yuMfdpxuPZwFH6MGqcIC3MWhhYRQujEnEBENLUdHxfu87VH2R5jUidQC6eVrHXrK8bS\nRbTcT9Yul1Mg0x0GlVqUUOur9+x1kRogMNPLPEGI3lATL54MR9C4fjktUBY3S1laUtvmiT61vNtB\n0TpZ1MMGvYXRE8fvFID5vxbxz5/5QSI0Bnwb8bUbg16c6/zyBvDcHeO5eHjJtWjiHBbJyzD+IA/j\nrsjQN/gLTOXvwnajdB41AwVIGa6vwoDA/I/OdY0c5jlTp8FOhvXNdp9Ur5rOTjcFKefYvBwwGPTs\nAuhAZhyBgkG0p/vVmQ7ZZLnUQ1G7tSAzQtOlIrWoQMPIVj7M4Lkhx6MPfdEcRP4/G7D9SUTTont5\n89rdi6NcxggEcbPTvdj+1GPti8czzqnx67jNTJlI02pRQ00ybqgG0Er6uugma39ZjkKBkOvXolUG\nulxWLDLIBtDtGqixuluPNUyK9jS/V++hIbvjIZxR1hBbw0xXrpUe9EQXLM9j/uck/Kcc2h9fd5wb\npqZ5JyXaUM4h/j2ZST7mUhjIZxgjhMty3l1LLR3I+J1eybxwjaTlupyuuAjCOTNg8Ro3ks8hfH2b\ncsyVra2xld5N4EFjb6bi2Fo1fpfL94aPNd0cENTa6OUSespn+rsc2ukSRsf/QXAa85/kchjLno1v\nAebjx3HtfoRrjjwihkgfIjIDFoFBDTfuOnQg69I3FY0M+Pgy0gMICqAOGGrkoKYJ3C7r4NIQbdsB\ngjo3B0YrTyBmF6YSFxjvDASa17lwKEs3nC5uglp5K4sOpyOANkw5FUc0LSu6uKPP6ZkGAOY/a2ZW\nP6KIAAVm8TfnH/EixlSAQWCGeZQALC7wDSoTD2nlY73XlHFALLXxVYfBc52tP243A4aaXi490Pp8\nXCNXbsf1yQbv+nYR8nY5QDDgENUZrw54FrLXwqEehNZQLSOKMvRXUHIpA1/TzJR17bpFxYtRp469\nSLwfw7U/TwEwevc7sAh6nN+HThFFTI9Lu8Y9xO7EMwA+M9U7g/GpSqcwpiQBLryL4SK4SA/cw0Jq\n4xOfa6CeeUOdW42ytC8XpTDPEcdY+BF5rk9eL8pH8Oca4ZmlpU4fp7uXFYgMjsGHh5rbujy3B82d\nEatxAnNP6PKuGmHFkm0T8WcXhYTX5UWYAVfNm8Q7cwAc9iuvEj9uUjLKyRbmvyaM1CB2DDYxv9U5\ndOY7+/Qfp7lN1p2vR8dSSUVNb1h6gKBFBmaAULv5KQAubrKK+gHIfG1Koupj/fjVIzvhTRblgCMD\nZft5AlzertLaKeA2gUVDyRaS64ONU/PJCJNr7WQTmhlILaXR+i3uQyOpi5gbKP9iMLw+l3cAfXl6\nj6bqJJ0AACAASURBVF9Axh+ixG8hOLXgPjj10DZZshQtpGUcPUvaEX7Z9cIc1ydc8drQeXOpBUcG\n3Bbrxy+nA9epcW79csCRwRbqOWNG4NUMJTNKl4q4CCL6VQb8MryeTOT1DGd2rbtH9kWA42tkvYLs\nC5DbxPyJSCeoLF9/LOC4jfgoxv9hiHLBPTCfE6lF8AP6xygcCiswq6FEP3zOpQjLLmeNNp24qDXq\n9PAXxLUs9LWJfK0GyLecXSsSzdLh/tauo2SG1Ed4jFLL/3UiHKHn6oYOGvbxpDNiZ/ljttWjn7l+\nja120+XquuuOBaZAGncqqrdk4pB/AcnE3yYWt32VELxC39XDafq2SS8XRfC1afjsOBtu24kDghoB\nvAxY19rMjrFeM/M5dHCpc8Z59PS3bKl9Fd4pYFlWNRdF1EK/EIfovFA1j9+JjryQtT0NI7m8Alot\nhFbQUMJoS8o6DxWgF78jOIZrP3deiCTiWlhf5Sci4nBz4Ag2bcPxNy58juvRMXI8TE/0kBmyRgiO\nZ6mtu5YoB6bXqetAr8Wl3ArwdQd7gGDgCECOFOJ7C0X5swOEzGvrMfZYNU/Pi4Bf3KaLUFy/uqA0\nTHbSc10xhgw2bHw6TpoaACO34MJb1jfK83iE8P0GTCpyJMEpSAABG0HGLzhA4jXD1w+p0yPZ2PaS\nlnE9Gfmb6a398bHWegjJUlk3notywJGBeqtafsXfncTFKmOuZXr0ybiI4DiAxX349kD7yQZ8asB1\ngO0LzUktZ3ULmI1QvfMM1/4m6do13jGVYYNg78seLeQocg+pc8+3k9fG0p3TY+pUasLgzYScemUV\nbl+jLAajnkhU9WBRHTjCyezC3Y278pwBixsMFx6z1KKIbEHpL/Y0/MsMR/XghavhrEP+zCgyrsTl\nhstMl0tveIw1JA/ZxPyhpLE7cC/m5GKMX6t9vQ42MMeu88uNWS1F28ky5nnWtaepq2vfAVBmzC4y\nnJlzmVOJuXLjoU6m5kDq4/SG6tm1rGUtrxtZkcigJ8/JVF2GqOGcTUm21jalpi96u27Whrs2R/C5\nOi5l2qmoV+F8X/tSPfgnz1nIHOU17ahtu8UcO7Is05ullTbGu7bncvOszWxdZLxLTS/dQeFjWVTs\nUp6ZvHN5/twTKc9lRcCgNhAKCNkWXQsUlJPgttwkOd00dHf8RK0t125GkrXq9fajwuG628feojLH\n6diL2L7vrylCLT/OQlklYTPj3+l9GJzH19rrMQXVWw1U2+A1p3zULKnj+uP6St5GOsNArOlGn5mv\nCBioZFuNytrzuxM3Oa6My9sz1AcdUwY72uvRLYt4nBfjSKJmcL0GE8x+1g6DZizauBWZx4MBge/L\n4DFw45HxPFy3FlG0PG+UaY2x6uiYeReNsN5uTXKfjjPIdA5hHkCNmp1AjathHfrkgMGAFXcLuZcR\nVjKLPViL1ImB18F3QOBIvUz3TMdeqYXjWd9ZtJARcgwKmSHq/RBOHwYE/s59h+hijTEOEGZ93L6+\nXhOfd2PsHvum4nRVYtnJsmtA+1PRbfYaiGR1QhQc2vqtUGSgoVwYNXuwWmisA+f+cKMmbDCK9j07\nAi1vpqLhXfSloltKqlMrv67lilFeHwga59TA4nxNZ/VWGSvv9NI+Wx4185Qu2mAdWHfuw/2CUq9f\nHURLV44ae6XGM9WEU8AMwHNZETBwJI+e44nIvIAeWxYQAI/2DrF7BpwXnVsUDAg1cYRXtmhrJJsD\nXNAxF/63xIEG651FBo4QqxlWLSyv6ZCRwPFnqjN69Rhc9Mvz7vp2c9viKbI02K3bGsi7dKV9XQcI\nBhqmu3OM9BlhVyPUapeXlXP5fi1crbXL7WtOrtLDIThvm5XltphE07LZotOIJFtMtUWWhfVu/Bz4\nZCma6peBV2YwnBr2pjcaObr10RsV1oy7NzXJ2ujVYbscIBi0+ADN2WIx840+boFxvd6HoXIfcZzb\n0YnnB4a4tpwXC9Gbdri86pGNUXubKCfpWNSjMchki1KNwY1Rlg44XdiQXf7rnqvouIpeI8gMhtvc\nDTfEbSiwOKDJ5pLLMjmtdTLPz/PSZ+bNm45KKb9cSjlXSvkMHfuxUsrTpZRPTa/30LkfKaU8UUp5\nvJTybfXW2eiyMPOKeTGRFhOmk3Y5qcseAdJOq3+XNszoXLxa3lINMPNGulB4l0PBMKuvoJDl6lvm\ns+owkzKhuwMG1XVZUV25rRrAcX0nbq3xddTmTtdZTRxZ69rbqWxg++PpmC9gHXqurW+WfgXAzwH4\nIB0bAPzsMAw/ywVLKQ8C+C4AD2K8h/XjpZS3DsPw2vZmXUioaMiGE985xOPJ2YA3NO7PPUgC5hiD\nVHzWhem2oVz4vowXd9K7VRSR0zKRhDOMlr7Zgso8q6vfWnYuz+8FFCULVUd9LFn0kRlrLW3jcq2o\nIUuN3HrRedQyzgkoL9VKKfs0XJBhGH6nlHLGnCrm2HsBfGgYhi0AZ0spTwJ4GMDv5l1ng8ghe48w\nSNS8MwNCtlhd24HE/K9NXN4ZVja8PGmZp3b9a/kM+LhOKwzl8q6OpgE9OW2NtNqUdxXntUOUe8lA\nXa+BvWKWljAhDGw3xi0p5/qN9lvGG+XUefHY1EyT9WTwyCLMPtkNZ/CDpZTvBvAIgP9+GIbzAO7B\nouE/jcX/1zbSEz45w3FekMPITFpcBbfBJOcGRiDgvywDtv+ngy6aGokV5WrA15Nf7pVwCM4G4NIj\nlWzOWuV6waoVUfQQqjVHwdGeA8YWj5M5nxrv4HiW2voNhxNrU8Ek9MzSvLrsFAz+AYD/Zfr8EwB+\nBsD3J2UHf/hXAVwF8BqA/2B6sUo8IXxPvGO8ncHUQreakWo7vJMRj+5izxevjKzMmOYayRj1lIRq\nEXyOdNI+XL1NbDeAWnjMemQhOadamuo47+v6yv6e3aURmXAEkO0QqX61Nlugpx56U44ryDpwdADi\ndM+iS/7+SQD/Vi/Cyo7AYBiG5+NzKeUXAfyT6eszAO6novdNx4x8LxYv0i2aeGfjVZafjUQHPsuf\nHUOrZeIc714EEDg9swXqJlHruJQj3jlq6MlfnWSAFDpoGuKMldM6F87y/CjnEror36Nth8ykLqQs\nt8d9R//83kpvnCNZxiwyp+PGzpF7WofXHJfTNIf71/Fg+drpFfJLpsx2jbullPKmYRi+MH39yxgf\nog8AHwHwa6WUn8WYHrwFwCfaXfPiUE8Qi1VDcWcwGurVHq7BfbO4BRGL3g24i05cOQUwLruM1HLo\nZQi7uB5eSD1huBuX+Ox2SeI95mKDXpnOAdYZwclkMOumZbOogPvh+hnJp9KKtBywuvC91bb24dIc\nnd9WipxLcyWWUj4E4JsA3FlKeQrAjwL45lLKQxhTgM8B+C8BYBiGx0opHwbw2KTlDwzDkKQJmVH1\nig6IEoPRfhx3niILFWuen+toauAWAffJ7z0gkBmoElrxno1hlqJEW1ou2z1RXfTYBrb/tkBzbE2r\nFMx7rwFo/4dBrW5ICwSdONCpieMd+Ht8dg5EwV5B0umfRZp1Kamt7qOUUoaRZ9Q0QVnSWvjoPAb/\nC7M+PivKu60lbgPYDgAq6vVUGIyysq0+tExP+egPRq/eCImN00UM2r8a+hb8vRlR17XNfWgbej0c\nLfZ4xBpn4+Z7WUPXflj0+nrWb7Z2Qty9NlrXpSPx/VswDIPbCdwxgbhHwqmBm8yWB3XkTRaCxgAp\n+PD5LB8OXZxh62LXxa3ektvTPvVaVI+aIfIxJqpaXpGvQXN8btMRWywRASjx6frgvrL5rUUJIVm4\nrB7X7QC41KAXCGrEtLbJUlvLzGNFW7zeNpEDSYtL6gO0AwaDEF7AerwlalQZB8Btam6r/enCr+XC\ntbqM0PzDGNbVAUKmB59zwkCmC4vfszayhRZ6aFheW3TMSWzJZ0g/2V+8LZML18azFU05cQbkxrQm\nPSDqIjj25uFI3G3sLrXQtvpThQMEgwzdHCFSuxC9BPZuWcShZbO2diq8ADj64b56FrALZaMtFRfZ\nsOhCVlDkPjQEdeLIuYyDUT6D9azN7WajXG08a9I6v0yKlRGctYiH67pyLp1wY1sDmj6egGVFIgOW\nzDu0Ls4NMqcMKmqccazlRTj6iEXOKUG8OLzPPOoyC1jFAZkuEk1R1GhrCzZbpFy/lTLxe4/U5jgb\nw2zOsmirJQowvcKhfO81ZxGcSwccMLtrynRu67SCYOA8obKjWlaFI4MrUtZ56FZ7qk+QNbwQNRJw\noJYBlLtW1k/TnyDYdDzc/w7EA2I4bI82uD+Xb7Oujm9xQJCF50r6ZeHrFXhduJ0tKZMZAOuUcQpa\nNvtea5/fM+5rN8IEYbTnxnt3/awgGADbDakVSnI9rZNtP3Fu5s61+tPQl/vPQCUWpZu0GoHKXgeY\nA4FubbLeAQwMJlnfXIdFo7IaGKjOLApOWlbD4lq6kx1jPd1xBQRNuzJQqfEsWkf1YhK111BdewHK\ntT8pzqTfxA8YDDJjArajrIamwOJgwZTVtmPyaoumZ0iyBZqBmIbzHFmojlo3RA2RfzYdcgxzz9pz\nwxXrHv1nXEPPbkq0pdcbc6npUxin/gQ8kxYvovqqZCDi5kKFgaQWodbWwTLCKWiLKHR9LE+arlhk\nkE1IeLSakUPOtwgnwHuDZRaRgpHm2rVc2pFqTg/1tNk4RBvhPeJ/Dlg47K/tGjgdstSkteg0PO8x\nlF7v19o9yVKXWoqh88LfW2vDRYa7Cd15Tmp6ujFchruYt3hA4hhblxtHWUeq9Ladna+Riyo1r+II\ntyyMVs5Ar8l5XTexWdga5R03MTNl+XsWGbjryHTQ9rPxccL9ZHNY24HRtnr607Y1Oozvbu25iKq2\n6+N06nE0bq7cLkNrXddlRSIDR8SoOHa1d4E5cRPvyte2trJJysBL2X5gMVSeIQcRboPbyTgRxymw\nOENyBJ9KFgVxv1q2R9QJ1PJsBYRe0THXfqNtFRdxsrG6P4tlENE2as5qJzsgNenfqVoRMAjJBskZ\nqMuvawvRtdHK71l08cRCzPJ8NepaqBzlNSxX0VtOVXjhcf7NnsRFXr3hpPP22rYj1WpcCl83A8Cy\nYW4Gotq/Ak0G/hkZyMJAoGOhbbuUxIX+/O70qMnuzHkFwMAx1XpOJ8gtrp7cTA05JtEtYjXijDV3\n/QYhVosMXF65JS/WQaOBbAwYBLKQO8sxFdjYUPmc/vu05tWZd83SwGiDdWh5fJ2XnhBe+80IXGcW\nbockXvrU7lY0CvOekdKtccjKHHoCMSS7uBrJxVJbAFGH26gx++y1uK1e1l3bU10yFr0FbjoOkWpE\nXf2RUNZGT5+84MIAaiQckmOODAUWwTkzjmz+tX0WHR8tp0DtIkvHTYQe/BBSlS0qm+nnxK3d3roh\nvfO6vedDIhqaArmxLRteZvV2wnrXIoIstdFXjTisCS/oWMSZgahO2bVw21zeleFj3J4zMjfWHFUE\ncGQ7AS0duA6LAkE2/xz5OR1cRNDie1rRQzbX2Xpy7SwfEYQcAjDgkNotDicOjTflu/P4WZsZS6sh\nbeZxFMAct5BNNgvvEixrFCw1zsFFTU5c6uGiIgdGWV7OhOBul2btGjLdapECsPi/jS5C1TFwfXP7\nmb6sj+OPamnIzmUFwIBZfZdvZoQbSxYVtIxDAaEVSrr+uHwLpXsNSPtRj5X9yo/rcVSh3EMrjdL+\nVVxKwwbg+IsQN79K7mXlFSgc0DCYLRMu13YYoi0FgVi7GRBofeVlsvYz/Vgch6HnlkktVgIMgDYg\nsNSIIpbM8DLJGH5ta9l8rsakOw+lC9kZdVwL8wTqqTRl4DIZ75LxLyqOBNTrURK0J+XhVCEDcgaY\nDAgzyQhfbd+NXVY23hUQep1Xdqw3yutpe79r7rkoIITU2GUmeXrz+xZBpWX1qTIw5VyfmT4t7xHv\nmRfk95onUWZegUB5F83xlwGBrB2nU0gtxWOHoGQr18301PZUelj6rG70yzq56E4BIYv6VK9ldeF6\nuzPnFQIDoB5K1ZhumDpKxsF85rrsPZXM08igFb1kequn1GsIXfhdjXQZj5HxIty/LmAu38pVNVJh\n0dRGDbB2HZkHzwytx7CcPvGZ29D1pPWWJY0dYLXWjgKrm4/d8QNOVgwMgNxjaXjIXk9zRfXEfCwW\nEnvIjNByk8wkXjYhrXZcCKqhfJQLI+oBgVreW+M5HMjUUgcWjmCyhe9I06irBDHrmwEptx3tOePl\nOdKIYybHkdTj/mHKOn1UWgS0llu2fdfW8rIiYODyVxY1ZBcucp6XkVdZqpFNPue9vRPnQkP+gQ/r\nru1qCM+G1JsftyKPKOMiAycZgNX0UXCNa+Lx1LnSspkeveLIPpcebGE7ENc4FQalLJrgOpBjGvmx\n88tS5OsjKwAGDi1rYWVmGFekbCsn16hA+4lyyzGy2wFL0xU+z31nCzHaiLrZzgn339Iv2uwFgSy9\nYsnC66wP13/N89XCbmeUtfFUcTpmzkH1VuPm45kD0Lpxbst8dvrslFysywqAQQ8JBCwakE6wSylq\nHrK1SNyk9C7yDdQNrofw0zLRpouKMm+k3ktTCAeYmmZlEZbWaY1rHHMP9YzvLmVyooBeE420VPbC\nqLI+XNsZkLK4tGXv+YGs5xUUHZAaEADbQSBDehcOcn3tExgfGJKJM/Rsh8EZck14DJjU7CHOat54\np1KLtEIX9+s9rsvkayvNUrDLeIFMHAfDEuvJzZdGbs44s+hDPb/WzeZG0xh2OjsBreXnekXAIIsO\nWnmpA4WMGQ85mhx3/47kvK4jwqJ8LVrgdxX2xC4aqZFd2s8yOx3cXpbDumvS6+aUq5U3c73op5eI\nzcDARTmtcVDOpmenw6WZy0S1PR4/AyA0ju9eVgQMgDogtHJlZrJrSKpeJ0Qnqpa7R1+OfXfhMczx\nKJ+Fy5mX6fEQCgi1saiNa83DZ/1ytOXApBaxZYaVEW7cpjoALqfAqm33kMxcXiOCVj0H4Nnugq5P\nl/btn6wQGACebc0QPiORnAeviaYWWf7b066GsppPczldrEo4qi6sQyY6ThpptFKM2nFu0/XbAo+a\nV8zqu34UcLJwP65DH/6yU55Ax87xVi7NVNEx6bnubO3tLZewYmDgJAuNw/CyUDLEDZiG+mwsGvIq\n+tfyXQcEGYFVI05b6YYro8ZVI7AyQ8nEhfKb6DPi2q5FFqm58a0tVbf96QhUrdMC1wxEeV575krb\nqxHU2XpvgVgrhWzLCoJBtnXUc6EcXvWSNw4IXMjr0DlLa3rIJRVtyy2SaH8Z78aGoilBbxSwie3R\nCo9T1p9+rkUHmrtnXEILQHVs+B6PbB6yNE516vH20Z5LZ5y4Mcj0zNKLFpD0yQqCwV5Iy4vU9noV\nBID5YnZPB+a60V4GBFvYbhhZLsz9R7ls0bYkC6F7+uX+9XuWxtS2VrVPx8vw59b5Gdqgz3xSJi6a\naEltnJzebkcopAYerTHcGzkkYFBb8OpBHDEW5dQonJdzf+PVeowZsGjUGRC4v4NvhZ3OQ7LH6o06\neqISNTQHXk7n0EXrajtZv06ch3f1WlFSFjE4LueolK0Bv64lxxfUdgZq1+64lmWBanlZUTDIwvxM\nephcnpgMaXkSY4IDCFy4CfleY8ajLV7UGVkZ16QerTcy0PIMGrGQGWCWCS2zyIfbUeCsGUMPIx9t\ncPtZutR7XUpWu/Hnd8cdxeeMVHWiEUgvIPbYww3HGYTs1VYKe73aBLHwQlCPHsbkcnknLipg8lMB\nRic9C3FdpMPCEVLGn6gu2k4ritHrZL1aQNATKmsd1c2Vyc63mPgN87k2vtqPRnAZ5xHta9TAa2HZ\nSGpvbGWFwaAlrfA4894qNc/UO8it3QI3yWEwR01515b73iKl9H8N3fYl66X6usiCJbu+XiDI2qwB\nN+sc0trazFIFresAMRMHFBpFaDs13iID++vDFwCHFgx61M6Mi0Mzzcect3W5YObZuF4gvS5a541n\n8J7BLRCWmvc4Kt+jT/19gEY7biG7rTXVIeMXVDcFglqEoBFNrW1uk8GkBZy7ldoOlYvIsvTGgfH1\nlUMKBjUvmknNeNlT1niBnkmqMddhDBlzzH20Fn7NywXguX4YdPg9G8/adbfAi9sI/fh7q6yG/QFm\nrRB6Gb5pL6TVVyvFALaTkpnsT1QAHFowAHY+2TV+INoNY67lldpGi7hyC1vbYABwDP5OJDOMlsdc\nhul3qYHb5Wm1y2Wz6EsjBddHBpTLRAc9ZTMP7z63wD1rl/XePyAADjUY1MRdliP/VNRj7pXU8tss\nGmgZqeqo+b22E9HPFnzbNTDcknIu4pgBuGza5V8yZvMS4rxjLyBrW1o324FwfEhrDnYCFLtJV/YX\nBELecF16Wcta1rLycgNFBj2XotFBzSPUcrgsjeB+XLvahtZhz11LEXi70O2NOy+UkXvsLV3K5ET7\njajgJTp3AovPgsjGuiXLlKvxHkA+Zy6c302a1trWbLXJ28LXj/u4QcCAL6O1PeTShWXy2ZbU9ulD\nXOjLhsj/meiujevV+nfH3eKKxVcDE1eeSb4rAC4CuIAxNcieG5G1r8RmkIXZ3Li7+2rXxoDgdFG9\neO7UMfSuk2zMa3Oqcv12GG4QMHCi+aHziLHo3FaUfq5FAtmCdxyA8gWZtHY3oi23uGv767qn7vqt\n1atJbbxrZWvf+fqdDj3e09XLdOWx671+d7MT67/sOGZ97C8gVDmDUsr9pZTfKqX8f6WUPyil/NB0\n/FQp5WOllH9XSvloKeUk1fmRUsoTpZTHSynftq/ad4kzDp0g3orj7b+asbYYb/YulzFfXJexHSii\nPW07dKvtxTMhyf1dwOhZ4939LsKJEo583Hm6OL4J4NT0OoF5tLMsYTaT1xV51drT/jboBSyOq7uD\nktO0LXNM50tf0Se36YC1pnOt7P5KKzLYAvDfDsPwaCnlNgCfLKV8DMBfB/CxYRh+upTywwDeB+B9\npZQHAXwXgAcB3Avg46WUtw7D8No+XoNRubaN5PiADXl3UYUL0TNWPcv/o+/L8M9W1Db1MwNQbZuT\ndXGiC1cX+Zacd3XceJwyfbqx65Ge7bRsnl0UV/Osyxij7j4sw1W1ytTkgNOEYRieA/Dc9PlSKeWz\nGI38OwB801TsAwB+GyMgvBfAh4Zh2AJwtpTyJICHAfzuvmh/TTIiLZOMIGRiznmOKJN5gWhbgUC3\n5oARENxWZxg9ewwHBK0wfxOLfS9jkG77k9vIIhUFYv2+TK68jNTupVhmKzJEeQbXbi1lUd6Kx0bb\n4bqZXB8SsZszKKWcAfAuAL8H4PQwDOemU+cAnJ4+34NFw38aI3isiLiFmgGCkyys5/NZWMniFiGn\nJ1vYbvAOCNgINqVcfOZdhyyyaBlpKypw3Ajfr9C6N0KPZURgFhWp0S0r7np6PHmLw8jOO8l2nFZs\nN2FKEf4xgP9mGIaLpZRr54ZhGEopQ6V6cu799Pmh6bUbaXmBELeVpqJpQ43EC9mSF3tovfnF9c19\n9hhPpqN6b/ZULg1hr82LrydFaJF+kO+1aEPFbZ1qzq9RTy06aKVMrWMsGZi7z7U2Wv3sBRA8Or3a\n0tS6lLKJEQj+z2EYfnM6fK6UcvcwDM+VUt4E4Pnp+DMA7qfq903HjHxvl4LLSY/RAts9iNvicm3U\nwjqOChxR57yZbl21mP5MlxCtdwzbr60V8rsUqHa85omXMUAGtYw3yTiVnRqNA6deaa0ZTSF1ze1l\nmlQTdbQfSEu2dhMKgF8C8NgwDP8rnfoIgO+ZPn8PgN+k43+1lHJTKeXNAN4C4BNL6b4nshtEzdh0\nbZ/Z7VruCWzfoXA7FtnWYaZj9uL+jsqL+3USuh2lcqp3SG2MM+DTZzq0pGWsDngzvkN3jPSl9fUV\nwjsdcU06/zrvPZzLwUsLDr8BwH8G4NOllE9Nx34EwE8B+HAp5fsBnAXwnQAwDMNjpZQPA3gM42j9\nwDAMtRRiH8WF2BlrHuImtVVGSTL18mqAjkxqed3dgFtt4Tly0eXkLNk2LacCjnyrRRCOu8lEd3nc\n55Dazoduj9aARz27bl9qxNi61tWUchC2OnIMv3Wde615NRe+uf3i7AaSDSnHbTIIZEaW9am69RiM\nGoBb2NpulOfoAVJOd0vcZ25Lx0F3YzJCsNdYMqCscR0s7hpcG1HWXWu0z+OWrZtVkW/BMAzFnVk2\nUTrEwjmdI91apCLXa5GAKj1RAOsZbWeGy+KM2oWjeo2921qhk4Ih75pwXxyGq0etpV6hR20rVPtm\nYXCJtvicHtM2ljXYFtj2yG7q7r28jsAA8CSPAwTAL/LaFlYt91dRL8JGr0aXGYfz2PGdw1719NpG\nRiCGKBAoCPACPobFa3NG664h+u3ZEXLeVuu0QCE+b2DxAbXZboRKNrbuWlsGvjpRw+sMDJzUAEFl\ns1HGGZ0LQdWLRbnsPoXwmM7Tq1HOsPgjId2lyBafEmnAdoCpAQHrxanAUcyJtmzsePxbgJDt1nDq\nkf1IyqVhPCaxFrI0MgOWVpSXyeoAAfC6BIPe+xGccOir4jxRDVhUdFdAPbyy2q5tXrTZNbrtL949\n4Ggldk0yA3Rthq7qnQMIsjy6di+Cnuf0I+McXP5eA0JgTnzWjLSHhzic8joEA2C7J2DRdABUlt9V\neodSwSQMg6OCXsNbBmyyNngbUolQl7YAeeoEKePSoZ1sMWq/EW0wWaqAkM2H4x16dp74nB7b6Tys\nlrxOwSAkW4g9jHZr6Hp3LzT0dnoo+67bmNFOTW+NatRrK8goQDm99DPX5+jEPQ6+9QtEpz/3mUUj\n0T/XyQjQZSXbqtxpe6slr3MwWEZq6UVG7oW0cspse47r9oBP6MjgwQbh2sjCab2pqrVFGu1zqO3A\nh/vQ62/tNABzYInoYBlxaUqvCWRAcGNEBcAaDJYUt4BaLD9LbbssO++4AhdNcLmjpp6TOK/byINv\nfgAABOBJREFUmfrsAAWCVrs6Tkfhx4O3ErWu6sj962e3+8LCkY5Le5YxAwbGGwcIgDUYLCG8CNjL\nthaEY7QzY8q8Ty3lcOGyRgVa1+leIwvVCJ24LUn9W7oQB4B8Xa30zf05LosDmNo9Ci5qqhGJLQA7\nnLIGg6VEc/VltpR6jIlJsVq7jmPQqKA3v2Xj1zxeuYoMlNzxaC/O90YVy2zbqsG2jLJ3N8SJm+8b\nAwRC1mCwK2ktEHdvALAY3sb3FgGYhcFcJwOC7E67OO4iAt2372XP40aeGeaPeIvjsWvBOmYAoykK\nk5IzelddYsyzrUtuU0XHa1mwOdyy/t+EbtmkF8sj8ATYMnvycdy9tG5234Hb2uOtQU4DeHuPny/4\nKWw3mjBg/uWjk4zbiLYvou85hjzOGVHJUYzjAbYAfJKuu3Z/QStiq92zcT2k71kEeyFrMNi1xGTt\nhdfoTS82zUtDec2Rw/gVFNgwH5d+M0NwhsrAw9/ZcLk/RxS6G6L4fLTr7stQecToFu0qCOzFVvJ+\nyfUDg3Wa0JTWHWmZ9CwwFhdd1G5+ibIuTeF6rZuYuO03YO79j8FHQ1kIzn06ktNtMTrJbiTifhRQ\nsmtbltztqXPjyhoMmtILBMoD9LZXW3y17UklG3UfXyOD7A5CvonpZmzP6aMv7ZfvJ9A+ua+j8tkB\nGtfRm6FUWgDgxJGdayBQOcDnGaxlLWs5CMmeZ3AgYLCWtaxl9WRNIK5lLWsBsAaDtaxlLZNcdzAo\npXz79D+MT0x/zXZopJRytpTy6VLKp0opn5iOpf87uUpSSvnlUsq5Uspn6Ngh+s/MuSTX8mOllKen\nuflUKeU9dG4lr2Xl/st0GIbr9gJwBMCTAM5gpHMfBfD266nDLvX/HIBTcuynAfyP0+cfBvBTB61n\novufxfiPWJ9p6Y7xvzIfnebozDRnbzjoa2hcy48C+O9M2ZW9FgB3A3ho+nwbgD8E8PaDmpfrHRk8\nDODJYRjODuP/Mf4jjP/PeJhEmdjvwPyfKT4A4C9dX3X6ZBiG3wHwp3I40/3af2YOw3AW46J7+Hro\n2SPJtQDb5wZY4WsZhuG5YRgenT5fAsD/ZXrd5+V6g8G9AJ6i709jpf6LsSkDxn+WfqSU8jemY9n/\nTh4Gqf1n5tNU7rDM0w+WUn6/lPJLFFofimtZ4r9M9+1arjcYHPZ9zG8YhuFdAN4D4L8upfxZPjmM\nsdyhvMYO3Vf9uv4BgDdj/C+xLwD4mUrZlboW/S9TPnc95+V6g8EzWPwvxvuxiHQrLcMwfGF6/xMA\nv4ExRDtXSrkbAOR/Jw+DZLrrPN2H9D8zV0OGYXh+mATAL2IePq/0tdT+y3Q6f93m5XqDwSMA3lJK\nOVNKuQnAd2H8f8aVl1LKsVLK8enzrQC+DcBnkP/v5GGQFf/PzH6ZjCbkL2OcG2CFr2Xl/sv0ABjU\n92BkTZ8E8CMHzeguofebMf9/6z8I3QGcAvBxAP8OwEcBnDxoXRP9PwTgWQCvYuRt/npNdwB/e5qj\nxwH8+YPWv3Et3wfggwA+DeD3J+M5verXAuAbAbw2ralPTa9vP6h5Wd+OvJa1rAXA+g7EtaxlLZOs\nwWAta1kLgDUYrGUta5lkDQZrWctaAKzBYC1rWcskazBYy1rWAmANBmtZy1omWYPBWtayFgDA/w/r\nzQAAcbED0gAAAABJRU5ErkJggg==\n",
"text/plain": [
"<matplotlib.figure.Figure at 0x112c01ad0>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"psf = py3shape.utils.getPSFExarray(psf_data, x, y, nx, ny, options.upsampling)\n",
"imshow(psf, interpolation='nearest')\n",
"title(\"PSF Image\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Run im3shape. \n",
"The returned fitted model is normalized to 1.0, so scale it up to match the input image. \n",
"You could also add a weight= option to py3shape.analyze with a weight image the same size as the data image - you should do this."
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"results, fit = py3shape.analyze(stamp, psf, options)\n",
"fit *= stamp.sum()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Show the fitted image and the residuals"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"<matplotlib.colorbar.Colorbar instance at 0x1139d10e0>"
]
},
"execution_count": 10,
"metadata": {},
"output_type": "execute_result"
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"/usr/local/lib/python2.7/site-packages/matplotlib/collections.py:590: FutureWarning: elementwise comparison failed; returning scalar instead, but in the future will perform elementwise comparison\n",
" if self._edgecolors == str('face'):\n"
]
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAToAAAEKCAYAAACCOJnxAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztnX+wXVWV5z9fAgGUIB3pCb+iRBLGyjRjENt0j7/iNGJw\nbMAei0BVA2p0rMoIjNVjQXSqBR1poRr81QU1NQRIUKMZf2BUkASHWFhTEMEEAuGnEpv3JAHlt/xK\nwpo/zr7k5L2977v7nXNu3jt3fapOvfPW3Wf/OPfedc/ea6+1ZGY4juO0mb32dAccx3GaxhWd4zit\nxxWd4zitxxWd4zitxxWd4zitxxWd4zitxxWdMyaSjpT0iqQxPy+SPiLpln70y3F6xRVdC5G0RdJL\nkl4/Qr4hKKw37Km+Oc6ewBVdOzHgt8DpHYGkY4D9w2uOM1C4omsv3wTOLP1/FrACEICk10laIemx\n8AT4OUmd1/aS9M+SHpf0G+A/lSsO1y6T9HtJQ5K+2Mu01nH2FP7hbC+3AgdKerOkKcAiCuUHhbL7\nBjANmAW8h0IpfjS8/l8olNs84G3Ah9n9SfAa4GXgKOBY4ATg4w2OxXEq4Yqu3VxLocDeB2wGhoO8\no/iWmtmfzOx3wKXAGeH1U4GvmNmwmT0JXMSuJ8EZwInAp83sBTN7HPgqcFqfxuQ42ey9pzvgNIZR\nKLpbKJ7aXp22AgcD+wC/K5X/V+DwcH4o8MiI1zq8MVz7aJjpQvGDWS7jOBMKV3Qtxsz+VdJvKZ7A\nPlZ66Q/AduBI4N4gewMwFM4fDf9Teq3DI8BLwOvN7JUGuu04teNT1/azGPiPZvZCSbYTWAV8SdIB\nkt4IfJpda3irgHMkHS7pz4DzOxea2aPAGuAySdOC4eIoSe/uy2gcZxy4oms5ZvZbM/t1WRSOs4E/\nUWxDuQX4FnB1KPO/gRuBO4Hbge+zuzHiTGAqxbrfE8D/AQ4ZUb/jTBjkgTcdx2k7/kTnOE7rcUXn\nOE7rGbeik7RQ0n2SHpR0Xp2dchzHqZNxrdGFnfb3A8dTbEL9FXC6md3b9ULHcZw9wHj30b0deMjM\ntgBI+g5wMrv2ZCHJrRyOswcxM41dKs54vr9V2mua8Sq6w9l95/wQMH90scvD358AHxxnU72yvXR+\nA8Ue2SbZUTq/EXh/l7Kp27xPQv5CQh5rb6yyZfbPKFu+n2spvMgg3ec66Nyn8vuXam97RJbbt04d\nVT4v47kf/fg+LKlcw//MKPs/KrfWLONVdD1q+5+Evw+E4+hxNuc4Tnc637H6aPInrd+MV9ENAzNL\n/89kl/tQic6v1k9wJec4TXI0u3/Hrq9cY5v8Q8c7ltuBOZKOBH5PEQnj9HTxfiu52X1u76iWt/em\nPrfX7/ev3+1Njh/9nIWOic64FJ2Z7ZD0KYrFoinAsu4W1zrf2NjazEjm1NheivKte/MYZXM/MrHy\n5bW48hcz9RbmtLkjIitPXMrjy/nIxOrtRqfNuZnXVaUfn5cyk0PR+dQVMLMbKFZxHcdpIT51dRyn\n9fgTneM4radNyqFNY3Ecp0b8iW7Sk7MJNbeOOmxVsbqfyexHTJ6zuTh3k3Md5GwCjslzP86x8ikD\nSs5nox24onMcp/W0aXuJh2lyHCfK3hnHSCRdJWmbpE0l2XRJayU9IGmNpIOC/EhJL0jaEI7LS9cc\nJ2lTiJL0tZJ8X0nfDfJbQzqAJK7oHMeJsk/GEeFqYOEI2fnAWjM7Gvg5pVwkFEFCjg1H2VH3CmCx\nmc2hcFLo1LkY+GOQfwW4uNtYXNE5jhOlyhOdmd0CPDlCfBKwPJwvB07p1r6kQ4FpZrY+iFaUrinX\n9X3gb7rV5YrOcZwoFZ/oYswws23hfBswo/TarDBtXSfpnUF2OLv70A+zK/fwqxGUzGwH8LSk6amG\nB8AYUYftqA4LXw4pC2sOORbWHFJL1LH7kepDHZbNHItpHe9VzucoNY7JZcds8hNuZlaKefd7YKaZ\nPSnprcB1kv5dne0NgKJzHGc8dFPL68ORyTZJh5jZ1jAtfQzAzF4GXg7nv5b0GwoH5GHgiNL1R7Dr\nCW+YIrH67yXtDbzOzJ5INexTV8dxouzf5XgP8JnS0SOrgbPC+VnAdQCSDg7pGZD0Jgol99uQLP0Z\nSfMlCTgD+FGkrg9TGDeS+BOd4zhRqky0Ja2k0IcHS3oE+Efgy8AqSYuBLcCpofi7gS9I2g68AnzS\nzJ4Kry0BrqHQr9eb2c+CfBlwraQHgT8Cp3XtT1MJrIv59+VjF8wmd4d6U2srTXoOxNboUutdTa37\npMZ3YEb5OtboUuNrao0uN7RUjImwRrekcs6IrRnlD6GdOSMmEXW47sRuU2pRPvZFSX3hU30bKz5c\nmWkZdaSI9aMpVzZozlCSog7llfqqxMaSam9yGSP26ffvQ4MMgKJzHGc87O2KznGctrPPlD3dg/pw\nRec4TpSsJ7oJTouG4jhOneyz757uQX1UUnSStlCYCHcC283s7XV0ynGcCUCLHoOqDsWABd12JE9M\nmtyGELPC5VhXU/1IWUFzAlOmqMPCmmNtThHrc1N9S9GUlb6u9vpouXVFtxsTdu+M4zgVaJGiq+oC\nZsBNkm6X9Ik6OuQ4zgRhSsYxwamqs99hZo9K+nNgraT7QhyqwE9KRY9msiTudZzJxwPhqJEWPdFV\nGkpwusXMHpf0Q+DtQEnRfbBK9Y7j9MzIB4nrq1fZIqvruKeukl4jaVo4fy1wArCp+1WO40waqoQY\nnmBU6eIM4IdF9BT2Br5lZmtq6dWr1JF+MGb5S1kDm8p7lGspmwhBHpu0BtaRKrJqkM5cUu2lghz0\nWsdECACQYBIosF4Z91DM7GFgXo19cRxnIjEJjAy90iKd7ThOrbRIO7RoKI7j1EqLtEOLhuI4Tq20\nSDtMkKHkBKbMXaTNcUOqw3WnSfefGE0ZLppcDM+5RznvX64xKef+57iRpcZXR4DTPgZ+a9H2kgmi\n6BzHmXC0SDu0aCiO49RKi6yunu7QcZw4FTcMSzpX0iZJd0s6N8imS1or6QFJayQdVCq/VNKDku6T\ndEJJflyo50FJXxvPUFzROY4Tp4Kik/QXwMeBvwTeAnxQ0lHA+cBaMzuaIhfr+aH8XGARMBdYCFwe\ncrkCXAEsNrM5wBxJC3OH4orOcZw41aKXvBm4zcxeNLOdwC+A/wycBCwPZZYDp4Tzk4GVZrbdzLYA\nDwHzJR0KTDOz9aHcitI1PdPwGl1Vq2JOPtU6rGIpYnXXEaQzp88p61xOP+rIRZtjIc9tr+rHMfV5\ny7FsPpuQ1/FVqcNC3keqDflu4EuSpgMvAh8AbgdmmNm2UGYbhSspwGHAraXrh4DDKW7EUEk+HORZ\nuDHCcZw4+43/UjO7T9LFwBrgT8BGipQL5TJWJLpvHld0juPE6WJ1XTdUHN0ws6uAqwAkfYniyWyb\npEPMbGuYlj4Wig8DM0uXHxHKD4fzsnw4Zxjgis5xnBRdtMOCI4ujw4W3jS4j6d+Y2WOS3gD8HfBX\nwCzgLODi8Pe6UHw18G1Jl1FMTecA68NT3zOS5gPrgTOAr9c4FMdxBprq2uF7kl5Psc62xMyelvRl\nYJWkxcAW4FQAM9ssaRWwmWLxeomZdaa1S4BrKBZbrzezn+V2RLvqqpdi7t3rlpc6snLluArlujdN\nyyhbRzas2OL59ETZOmL25dRbx33udwy9nH6kPospg0ZMnnOPcjPE9fp+n4uZjTtxlSSzz2aUv4hK\n7TWNP9E5jhOnRdqhRUNxHKdWWqQdWjQUx3FqxaOXOI7TelqkHcZ0AZN0laRtkjaVZEnHXMdxWsKA\nZQG7GvgGhY9Zh45j7iWSzgv/nz/60l5dnHJcdJ7JKAt5bmQ55FiKc93WYlbenGCOKVLWwKY+qXVY\nV5vM4FVHe7G6c63sOcTuXUPuYoMUpsnMbgGeHCFOOeY6jtMWBuyJLkbKMddxnLYwCRRYr1QeSnfH\n3BtL50cBs6s25zhOlAcpIhvVSIumruNVdCnH3BG8f7z9chwniznh6JDtJTWaCtFLJhrjDby5msIh\nF3Z3zHUcpy0M0hqdpJXAe4CDJT0C/CMQdcwdR/XjINcfMqcPqbpj1tE6fExT1uY6/GVzgnrmWPKa\nDB4Ze6/qSLuYM5ZU2Rxrf5N97iODNHU1s9MTLx1fc18cx5lITIIntV5p0VAcx6mVFmmHFg3FcZxa\nGaSpq+M4A0qLrK4NK7peF3tzFovryOCVoqlF5AMTZVPGiKaCaeZk38q9zzFjSaqO1LhfE5E9nygb\ne19zjEmQF3gzRb8ze+X2rwL+ROc4TutpkXZo0VAcx6mVFmmHFg3FcZxaaZF2GK9nhOM4bWdKxjEC\nSf9W0obS8bSkcyVdIGmoJD+xdM1SSQ9Kuk/SCSX5cZI2hdd6zbi1Gy3S2Y7j1EoF7WBm9wPHAkja\niyLp9A+AjwGXmdll5fKS5gKLgLkUeV1vkjQnpDy8AlhsZuslXS9pYW7Kw4YVXa9WrZygknVYnXKD\nWOZY1mIWxdzbnGNtznEvy7G6pqjDip0TaLUOS3gdboAThdhYGgr0WV/OiOOBh8zsEUkCYmkRTwZW\nmtl2YIukh4D5kn4HTDOz9aHcCor4l1mKzqeujuPEqc+p/zRgZTg34GxJd0paVkrDcBgwVLpmiOLJ\nbqR8OMizh+I4jjOaLtph3f8rjrGQNBX4W+C8ILoC+EI4/yJwKbC4Qi97whWd4zhxumiHBe8ujg4X\nXpYseiJwh5k9DmBmr8aulHQl8OPw7zAws3TdERRPcsPhvCwf7qn/JXzq6jhOFJvS+9GF09k1bSUE\n6u3wIaCTXXA1cJqkqZJmUUQRXW9mW4FnJM0P63tnMI74lw0/0TXhrpLqcmpBtqnF5ZwF9dzF4ljd\nMfcoyHOzyulzjttUbt2pOnLcunLqTRG7dznZ3bqVj5ET/y6nH818xndW1A6SXkthiPhESXyxpHkU\na3UPA58EMLPNklYBmykGvyRYXAGWANdQfMiuz7W4AmhXXfVS5JG4tMfSdejbHEXXpBU0h5RymIyK\nLucL31QdqXpzfmiaHHcd9faqnP8BM4tZN3tCkr34p97L7/daKrXXNL5G5zhOlJf2nZpR+uXG+lEH\nrugcx4myc0p7wpe4onMcJ8rOFsVpGtPqKukqSdskbSrJRvqrLWy2m47j9JsdTOn5mOj08kR3NfAN\nCteLDkbEX2181eeSWqTNcTdqyoIG8cXlXFeomOFhWmY/cowRTRmDUvczlVGrqivTnsic1d5J0c4W\nja2XLGC3SDoy8tKEtbA4jlOdgZq6diHmr+Y4TkvYyZSej4nOeJ9Ne/RXu6F0Pptis7PjOPXzEPCb\nWmt8iZztJRObcSm6Lv5qIzgxLnYcp2Zmh6PDmso1DtQaXQxJh5rZo+Hfsr+a4zgtYTJMSXtlTEUn\naSXwHuBgSY8AnwcWxPzVRhOzuuVYJXPIcYVKWedy/DhzLHy5QTNzfn9SZadHZCmr6+sz2ksRux8p\n62qqH0/0WC/As2P2aBc5FvJcy3SOtTknrWTOZyDnXvTOQCk6Mzs9Ir6qgb44jjOBmAz743qlPZNw\nx3FqZeDX6BzHaT8DNXV1HGcweXnQt5f0Tq8uPakF4JyMYTmx1vpNrrElJwDlgQl5JH/IAQlnltkR\nWWoL+IsJ+dZI/7akjBw59yNl0OjV0NWtvTqMPlVj3eUaqmJjbOZr7Gt0juO0Hl+jcxyn9fganeM4\nrccVneM4radNa3Se7tBxnCgvs2/PRwxJB0n6nqR7JW0OKQunS1or6QFJa8qRjyQtlfSgpPsknVCS\nHydpU3jta+MZS8NPdL262KQsTDnWuZyyucE7cwJy1pEGMVZHKvDmjLj4oIiF9a8SVbwzIvuLRNnn\nEvK7I7J1ibK3pizFOZm9cqyuORnimsxyVkdKyKrBSXunhqnr1yjSE35Y0t7Aa4HPAWvN7BJJ5wHn\nA+dLmgssAuZSbBm4SdKckPLwCmCxma2XdL2khbkpD/2JznGcKFVCqUt6HfAuM7sKwMx2mNnTwEnA\n8lBsOXBKOD8ZWGlm281sC0Xcqfkh4fU0M1sfyq0oXdMzvkbnOE6UittLZgGPS7oaeAtwB/DfgBlm\nti2U2cauKclhwK2l64conuy2h/MOw0Q3iXbHFZ3jOFG6TV3vX7eV+9dtS75OoVveCnzKzH4l6asU\n09RXMTMrEt03jys6x3GidFN0sxcczuwFux6sfnLhXSOLDAFDZvar8P/3gKXAVkmHmNnWMC3tBPEd\nBmaWrj8i1DEczsvy4dyxNKzoejU85Lju1BG/qw5Sfc6JOZZanI6NO1U2sbA/LyL7SLzo2xf9YpTs\nBG6Mlt3CrKj8uj+NXjZ5bsefxxsciosZyjFUxcgxGNRFU21WdS2rThVjRFBkj0g62sweAI4H7gnH\nWcDF4e914ZLVwLclXUYxNZ0DrA9Pfc9Img+sB84Avp7bH3+icxwnykuJbSMZnA18S9JUioQWHwWm\nAKskLQa2AKcCmNlmSauAzRSae0mwuAIsAa6h+LW/PtfiCq7oHMdJUHV7iZndCfxl5KXjE+UvAi6K\nyO8AjqnSF1d0juNEcRcwx3Faz8C4gEmaKelmSfdIulvSOUGedONwHKcd7GTvno+Jzlg93A582sw2\nSjoAuEPSWopFxVFuHKMv79VilrIaxeQpV6icunPfmJxAmDkWsFQdMetjIohl6icmYnWds+jOaNEf\n8HejZIefF8vIRZEPLsKFH3holOyCd14cL3xTXMxQjoU1JwtbTh25FsycDHExS2puUNZY+YlndZ1o\ndP3Gm9lWYGs4f07SvRSm35PY9ZFfTuHVGFF0juNMVgZG0ZWRdCRwLHAbaTcOx3FawkuDljMiTFu/\nD5xrZs9KuyJjdHfjKG86PYp4cgLHcarzEMVWtfqYDGtvvTLmSCTtQ6HkrjWzzi7mbQk3jhG8v65+\nOo7Tldns/iCxpnKNbZq6jmV1FbAM2GxmXy29tJrCfQN2d+NwHKcl7GRKz8dEZ6wnuncAfw/cJWlD\nkC0FvkzEjWM0MetTjpUp1r2UhSknuGKKVN0546gjsGisjoSva+odPGK06Hh+Hi16+OdGW1gvuCRe\n7ef+V1z+/qdG+8Ze8OaE1XW/uDiPHJ/pHAttbuDN2BuQk76zDuuxpzsci7Gsrr8k/dQXdeNwHKcd\nDNQaneM4g8lkmJL2iis6x3GivDxo20scxxk8BmaNrjq9Bs7McWHJ7XJssTfXcJETADSHVB0ZC9QH\n9F70+dQi+ZzRooR1iX1OiMs3M3e0cLRX2BjkZAFrKjBlrktWVfer3PaaqmM0vkbnOE7r8TU6x3Fa\njys6x3Faj6/ROY7TenyNznGc1uPbSyrRlBtMyvJUR90pl54YObc0x2KXsDI+lRj3ltGi7z69KFr0\nbR+5Y5Rs0Ue+Gy37A94VlX+eC0cLUwE2k9bY2BhT1tU63Lpi5AZUrSN4Z6/1pmjG6lrH1FXSFOB2\nihyvfyvpAuDjwOOhyGfN7IZQdinwMWAncI6ZrQny4yiygO1HkQXs3Nx++BOd4zhRapq6nkuRwrAT\nGtyAy8zssnIhSXOBRcBciuC+N0maE1IeXgEsNrP1kq6XtDA35WHX6CWO4wwuVaOXSDoC+ABwJdAJ\nYqnSeZmTgZVmtt3MtlA8988PYeCmmdn6UG4FMDpb+hi4onMcJ0oNYZq+AnwGeKUkM+BsSXdKWlZK\nrHUYMFQqN0TxZDdSPhzkWfjU1XGcKN320T277tc8t+7XydclfRB4zMw2SFpQeukK4Avh/IvApcDi\nyp0dg4YVXT/dvXLjiFVtL1U2x10sp71n4uKnDozLbx8tevFfpkeLnv3hK0fJPj87YlwAnrgv8WP6\nnYgsZYz4Q0IeHWNOhrhcw1OO8aMOcr4PTRrXeuMl9k2+NnXBXzN9wV+/+v+2C5eNLPIfgJMkfYDC\niHCgpBVmdmangKQrgR+Hf4eBmaXrj6B4khtm9+iKRwRZFj51dRwnSpWpq5l91sxmmtks4DTg/5rZ\nmWHNrcOHgE3hfDVwmqSpkmZReGCvD5kIn5E0P0Q8P4NxRDT3qavjOFFqdAETxdocwCWS3hL+fxj4\nJICZbZa0isJCuwNYEiyuAEsotpfsT7G9JMviCq7oHMdJUJcLmJmto8j9jJmd0aXcRcBFEfkdwDFV\n+jBWcpyZkm6WdI+kuyWdE+QXSBqStCEcC6t0wnGcicdO9u75mOiM1cPtwKfNbGPI7XqHpLUkNv05\njtMeBiZ6SVgI3BrOn5N0L7v2sMQ2/Y0gZsGKWZNSLlYxeY4rTqq9XJeZqtbYXLeiHGtgwhp7X8Qa\nm2ruvtGiJw5IWFeH4mI2ZpTl4YQ8lh44Mb6ovA6Laa77Vs5no9fvA+S5HTZjKR4YRVdG0pHAscCt\nFGkQz5Z0JsVGhn8ws6ea6KDjOHuGl15uj1N/T9tLwrT1e8C5ZvYcxaa/WcA84FGKTX+O47SInTv2\n7vmY6IzZQ0n7AN8Hvmlm1wGY2WOl18ub/kawtnT+JuCoCl11HCfNQ8Bvaq1x544BmbqGDXrLgM1m\n9tWS/FAzezT8W970N4L31dNLx3HGYHY4OqypXOPAKDqKtbi/B+6StCHIPgucLmkeIzb9OY7THnZs\nHxBFZ2a/JL6Od0Nv1edYWKuSYy1LlU3djj1vAUvX+8e4+KlInzcmLHxbI7LULUqZnF6MWb23JArn\nuCrWEcQyRY71vY4AoLGxpD5bOf7RzayRvbJz4q+99Up7RuI4Tr0M0NTVcZxB5cX2qIf2jMRxnHpp\nctWgz7iicxwnjiu6XokttMZkdbhv1RF4sw5DSay9VB9ek1Fv7mJ4pM0dM+JFh3Lc1p5IyGMuWdsS\nZVNjifUjZ7E/9+Mce1/qyAKWqiMmz3EDTNFMFjBXdI7jtJ/+BTNuHFd0juPE2bmnO1Afrugcx4nj\nU1fHcVrPi3u6A/Xhis5xnDj+RFeFnBXOmOUp101rTwRjHElqzM8n5FXTREJ83CmLadXgkSl5bpDU\n1P3otY5c176cOnI+t3Ws4ufsRGiICh97SfsBvwD2BaYCPzKzpZKmA98F3kjhI3hqJ5alpKXAxyhW\nB88xszVBfhxFcpz9KJLjnJvbH0936DhOnB0ZxwjM7EXgvWY2D/j3wHslvRM4H1hrZkcDPw//I2ku\nsAiYCywELg/Rk6CIf7nYzOYAc8aTo8YVneM4cbZnHBHMrPOIPhWYAjwJnAQsD/LlwCnh/GRgpZlt\nN7MtFAH25oc8sNPMbH0ot6J0Tc+4onMcJ87OjCOCpL0kbaTYOX6zmd0DzDCzzk7ybUBnF/th7J5l\nZIgiP81I+TC78tb0jBsjHMeJU3Fp2sxeAeZJeh1wo6T3jnjdJFn86nppWNH1aghIGRLqcP+JLd7m\nLhbnuPnEytbhmpS70N5UXLzUvavDOJDT51jdue9rTh25roe9tjeBnzW6bS+5fx08sK6naszsaUk/\nBY4Dtkk6xMy2hmlpJy3DMDCzdNkRFE9yw+G8LM8JaAj41NVxnBTdjA9HLYATL9h1jEDSwZIOCuf7\nU+RV2ACsBs4Kxc4Crgvnq4HTJE2VNAuYA6wPKVefkTQ/GCfOKF3TMxP458RxnD1KtanrocBySXtR\nPFBda2Y/DykZVklaTNheAmBmmyWtAjaHlpeYWWdau4Rie8n+FNtLfpbbGVd0juPEqaDozGwT8NaI\n/Ang+MQ1FwEXReR3AMeMvzdjTF0l7SfpNkkbJW2W9E9BPl3SWkkPSFrTeUR1HKdFVNxeMpHoquhy\nN/05jtMiKm4vmUiMOXXtsunvPUG+HFhHVNn1OjNOWdt6DdwJ6Z+VWB+mJcrW4dyX46JTRzDHnDqa\ndCvKcQFLfS5i5Zt04cu5z3W0lxPgNEWsjoYs7C1y6h/T6pq56c9xnLZQwQVsotHLE92E2fTnOE4f\nmQRrb73Ss9W1x01/IyjnuZ5NsTXGcZz6eZDCPbRGJsHaW690VXSSDgZ2mNlTpU1/F7Jr09/F7L7p\nbwQn1tlXx3GSzGH3B4nsrWajmQRT0l4Z64kua9Of4zgtYlAU3Xg2/e1Or36mqcWAnPSDqXclp44c\nf9kmfUxzFkfq8MHMsT7WMb5UHXX4r+aQY/GsI51mzueojjSIFRnENTrHcQaMl/Z0B+rDFZ3jOHEG\nZerqOM4A41NXx3Faz6BsL6lOrwunTS2oQ3yIOUYHiBs0cgwoOZmzIC8LWM4ieR3ZvupYlE/xTESW\nek9y7lFq3HVkW4tRhxFtcmcBm2j4E53jOHFc0TmO03p8jc5xnNbj20scx2k9PnV1HKf1+NS1V3p1\nLcrpRh0p7Zp0o8mx5KX6kQoMmlNHrM0ca2AdKSFzLYp19CNGbtDSnDpi97QO97k63NMq0qLtJZ7u\n0HGcOBUCb0q6StI2SZtKsgskDUnaEI4TS68tlfSgpPsknVCSHydpU3jta+Mdiis6x3HiVIswfDWw\ncITMgMvM7Nhw3AAgaS6wCJgbrrk85HAFuAJYbGZzgDmSRtbZE67oHMeJUyELmJndQpFfZiSKyE4G\nVprZdjPbQhFBdH4I6jvNzNaHciuAU8YzFFd0juPEeSnj6J2zJd0paVkpTephwFCpzBBweEQ+HOTZ\nNGyMqOquErs+x20qRa4xoik7e8o4kJNRK4dct64YOQaGPsZOA/JcvSDvntZhMKvj69bHPR/1N3UF\n8IVw/kXgUmBx7a1E8O0ljuPE6fb79/I62L4uqzozezW3jKQrgR+Hf4eBmaWiR1A8yQ2H87J8OKvR\ngCs6x3HidNteMmVBcXR44cIxq5N0qJk9Gv79ENCxyK4Gvi3pMoqp6Rxgfcgw+Iyk+cB64Azg63mD\nKHBF5zhOnApTV0krKZLcHyzpEeDzwAJJ8yisrw8DnwQws82SVgGbQ6tLzKyTQnUJcA3FOs/1Zjau\nrD/aVV+0s/sBvwD2BaYCPzKzpZIuAD4OPB6KLh3ZgSLX66U9diMnbFLuGl2veSsgvWbWx9A4QHMZ\n5HPaSzFR1uhyNp6n3tem1uiapNe8GudiZjELZ09IMvbPSNf8giq11zRjJcd5UdJ7zex5SXsDv5T0\nTnbth7kQ/mcdAAAGGElEQVSsL710HKf/DJILmJk9H06nAlPYtTemB+1dNbBhHS4zOb/ETVm06uhb\nHZ+6Ou59zlNTrutbrI6ce5S6zzlPw3V8BuoIJNvvWUSEFjn1j7mPTtJekjYC24Cbzeye8FJsP4zj\nOM6EY0xFZ2avmNk8CtPuuyUtoNgPMwuYBzxK74txjuM4fafneYGZPS3pp8DbzGxdRz5iP8wIbiid\nz6awGjuOUz8PhMOJ0VXRSToY2GFmT0naH3gfcKGkQ8xsayhW3g8zghPjYsdxaubocHS4voY622ON\nGOuJ7lBguaS9KKa515rZzyWtiO2HcRynTbTHGjHW9pJNwFsj8jMbarJC2Tr6UEcgxn7/CtaxF6yO\ncTe1Jy3HYp1rKW7Swtpr2VSfJ0DgzQF6onMcZ2Dpd1CG5nBF5zhOAn+icxyn9QzIGp3jOIOMP9HV\nTL/dXZpcAM5Z16jDmT5nob2ONZc6xt1vd71nM9qrIwBAzucoNebUfe7nU5Y/0TmO03r8ic5xnNbj\nVlfHcVqPT10dx2k9PnV1HKf1+BNdj0yA4IGNWUFz6k1Z/VLt5YQJz7HONWnJi9V9YKJsTtj7HKty\nrntaU5bpFP12k6uKP9E5jtN62vNEN2bgTcdxBpXtGcdoJC2UdJ+kByWd15cuJ/AnOsdxEox/Gi9p\nCvAvwPEUSad/JWm1md1bU+eycEXnOE6CSmt0bwceMrMtAJK+A5wMuKJrhpy8riliv2w5v3Z13Obc\njFoTgdR9zjWs9FpH6h6l3Lpi7eWuS1U1BEzkBf9Ka3SHA4+U/h8C5lfqTgUGQNE5jjM+uinhMXNU\nZGS/bp4+GSP6nbSj3+39ts/tPdTn9vo9vgdb3t5kSWKzo8vxJmBh6RjFMDCz9P9Miqe6PYIrulp4\nuM/t/abP7fV7fP1W5P1ub7IoukpW19uBOZKOlDQVWASs7kOno/jU1XGcBONfozOzHZI+BdwITAGW\n7SmLK7iicxwnSTVDl5ndwO7JnfcYMmtmzVDShFqMdJxBw8w03mvH8/2t0l7TNKboHMdxJgruAuY4\nTutxRec4TutpXNH127FX0hZJd0naIGl9A/VfJWmbpE0l2XRJayU9IGmNpIMabu8CSUNhjBskRTcy\njaOtmZJulnSPpLslnRPkjYyvS3tNjW8/SbdJ2ihps6R/CvKmxpdqr5HxOWkaXaMLjr33U3LsBU5v\n0sws6WHgODN7oqH63wU8B6wws2OC7BLgD2Z2SVDmf2Zm5zfY3ueBZ83ssjraKLV1CHCImW2UdABw\nB3AK8FEaGF+X9k6lgfGFNl9jZs9L2hv4JfDfgZNo7v2Ltfc3NDQ+J07TT3SvOvaa2Xag49jbNI1Z\nf8zsFuDJEeKTgOXhfDnFl7XJ9qCBMZrZVjPbGM6fo3DAPpyGxtelPWjoPTSz58PpVIr9XU/S7PsX\naw8a/Iw6o2la0cUcew9PlK0LA26SdLukTzTcVocZZrYtnG8DZvShzbMl3SlpWZ1T5Q6SjgSOBW6j\nD+MrtXdrEDUyPkl7SdpIMY6bzeweGhxfoj1o+P1zdqdpRbcn9q68w8yOBU4E/muY+vUNK9YCmh73\nFcAsYB7wKHBpnZWHaeT3gXPNbLc48E2ML7T3vdDeczQ4PjN7xczmAUcA75b03hGv1zq+SHsLaPj9\nc0bTtKLru2OvmT0a/j4O/JBi+tw028J6E5IOBR5rsjEze8wCwJXUOEZJ+1AouWvN7Logbmx8pfa+\n2WmvyfF1MLOngZ8Cx9GH96/U3tv6MT5nd5pWdH117JX0GknTwvlrgROATd2vqoXVwFnh/Czgui5l\nKxO+jB0+RE1jlCRgGbDZzL5aeqmR8aXaa3B8B3emiZL2B94HbKC58UXb6yjVQG3jc7pgZo0eFFPI\n+ylCRCxtuK1ZwMZw3N1Ee8BK4PfAyxTrjx8FpgM3UYSlWAMc1GB7HwNWAHcBd1J8KWfU1NY7gVfC\n/dsQjoVNjS/R3okNju8Y4NehvbuAzwR5U+NLtdfI+PxIH+4C5jhO63HPCMdxWo8rOsdxWo8rOsdx\nWo8rOsdxWo8rOsdxWo8rOsdxWo8rOsdxWo8rOsdxWs//Bw4tW6u9Pr0vAAAAAElFTkSuQmCC\n",
"text/plain": [
"<matplotlib.figure.Figure at 0x113412550>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAATwAAAEKCAYAAACPJum2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztnXmYVNXRh99yVEQkEDCgLBEluKAoCBEXlBEVUcPmGggI\nikrEhSiuaHTARAT3uBA/NxCEiAqKERW3UdDIJiAybKIkLAICioCKLPX90XeYgVN3mKG7B6en3ufp\nh5lfV99zT3dTc++pU1WiqjiO45QH9tjdJ+A4jlNauMNzHKfc4A7PcZxygzs8x3HKDe7wHMcpN7jD\ncxyn3OAOr5wiIn8SkbeKeD5XRHqkYJxsEVmc7HEcJxW4wysjiMgiEflBRNaJyHIRGSYiv9rV46nq\n86p6ZlEm0cNxMgZ3eGUHBf6gqpWBY4BGwO2795Qcp2zhDq8MoqorgPHAkQAicryIfCwi34rIDBFp\nmW8rIt1FZKGIfC8iX4pI50L6hEJ2Z4jIXBH5TkQeAaTQczkiMqzQ7/VEZKuI7BH9fomI5EVjLBSR\nK+LOXURuFpElke1cEWmVwrfGcYrEHV7ZQgBEpA7QBpgkIrWBfwP9VfXXwA3AyyJSXUQqAQ8DbVT1\nV8AJwIzgoCL7Ay8DfYHqwELgpEImO7u1XQGcE41xCfCgiDQxxjkMuApoFtm2BhYVc+6OkzTu8MoO\nArwiIt8D/yPhlP4OdAHGqeqbAKr6DjAVOIeEo9oKNBKRiqq6QlXzjGOfDXyuqqNVdYuqPgQs32Hs\nWFR1nKp+Ff38IYmrz5MN0y1ABeBIEdlLVf+nql8W9w1wnGRxh1d2UKB9dGWUDbQCmgIHARdEt7Pf\nisi3JK7ODlDVH4CLgD8Dy0Tk39FV1o7UApbsoBU7sioiZ4nIJyKyOhr/bBJXittPQPUL4C9ADrBC\nREaKyIHFHcdxksUdXhkkuop6BBhI4mpvmKr+utCjsqoOimzHq2pr4ABgLvCkcchlQN38X0RECv8O\nrAf2LfT7AYVsK5C4HR4E1Ihuq8cRc1WoqiNV9WQSjlqjOThOqeAOr+zyEHAcMBFoKyKtRSRLRPaJ\n9r7VFpEaItI+WsvbBGwgcVu5I+NI3GZ2FJE9gWsp5NRIrPudIiJ1RaQKcGuh5/aOHquArSJyFom1\nuQAROVREWkVOciPwU8z5OE5acIdXRlHVVcBQ4HqgHYmAw0oSV3x9SFxh7QFcBywFVpNYV7sy/xDR\nI/9YFwD3kHBcvyPhSPPHegd4AfgMmAK8Vui160g4yFHAGqAT8OqOpxv9WwEYAHwDfA3sz/bO03HS\ningBUMdxygt+hec4TrnBHZ7jOOWGXXZ4ItIm2im/QERuTuVJOY7jpINdWsMTkSxgHnA6iQXxKUAn\nVZ2T2tNzHMdJHXvu4uuOA75Q1UUAIvIvoD2wzeGJiEdDHGc3oqpFZsgUxa78/01mvNJiVx1ebbbf\nib8EaL6jURf9PwBm5ozlmJx2APyJ54ODHco8c5BFHBxot/M30/bjz0/b9nPO45DTK/Fzz6MeCmyf\nuO0v5jE4zZY/aHVcoC0utC93dE4e5+Y0BKAq3wa2p2yZEGgA92bdYOqr2T/QBs++vuCXx3LgqhwA\njjjyU/MYeVObBtp7zU4wbcfQIdBuKbQf+P6cH+iTk9h3XPvzNeYxGBJK5983LBSBFdQ09Qk3JLbv\n5XwMOSdG4uH2cMMvOy/Q8mho2g54tb+p68eJ/585EyGnRSRebY/3cF27HkJ1Vgdal/tfNm379vkr\nABNyPuDknER9h+ZMMm07XGmXKlw4OExMmRT+16Oz7LgzqOTY/9NsykrZnl11eMXy/jNzxgKwInce\ny3PncUC2ldXkOE6y5OWuIi93VUqPuVdKj/bLYFcd3lK2Tz2qS5iLue2qbmbOWHd2jpNGGmbvT8Ps\ngjuD0f3su6aSsKvO4ZfMrs5pKtBAROqRyMO8iMQOe5Oapezssn9fqsNxRPZvSnfA32eX6nAnZJfu\n3/rsuju3Sel4vy3d8X6bfVDpDriLVNzdJ5AGdjnTIsqZfAjIAp5W1QE7PK8DNVwQmcApgfbaRxfa\ng0w1tONjTqiCLd/cOCfQqmNf+r9GO1OfOPCMQFtyc1AMBIBVYZEQarLStL3czOOHm7kn0H7Nd6bt\nS5xv6tZ61ndUNW3fmhqu4Q1oZq9z3jo3XBMF7Pd/mW0q+9vfuSsPeyDQFlLftL2MpwLtH1xbbFuA\nRdQLtDF0NG07G2vPAKcQrs+ecHlQchCAHk8+Gmjn85Jpe+KWj039V3M2haJx2SJHJB+0sL+dNpeT\n2UELVPUN4I0UnovjOL8g/JbWcZxygwctHMcpN2Sic8jEOTmOkwIy8QovbeWhRETVWP9dOrpaoD1P\nZ/MYrxgLyB+/Z+8O/m8rO1J6UNdvwnPLtuc8vofVhgHqGX1mDv2rXQF9wl3hht8Ww+3NwWyw5Td7\ntgy0q3jMtF046ChT73vTXwOtHWNN26e4LNCu4R+m7TSamfqll48MtGeetAP3lw4PbQGqdwp2NvFo\n1lWmrRWUOQl7ob/NbR+YOrUMbYhtOnvKIaZ+5G1hSw45Nub/VLMw4LD6ILu18GLsUHVdo/L+SmMj\n9xHy36SDFva3xaYdGR60cBwns8nEbSnu8BzHMclE55CJc3IcJwVk4hqeOzzHcUwy0Tlk4pwcx0kB\nmXiFl9YobQsdH+jzCPNqZ3OkeYxWvBdo9fnCtLVSfABe4KJA+z8uN20P3rLI1O/JuiXQevJP09Yq\nabUvP5i2xzWZZeq8FH4mv6lvR4VXfhWTl2mleg2yTc986JVAe2fF6abtgJrhewFwOu8GWm8eNm2t\nqDfA8E7h53LcSDvCOooLAm1OTHmox7AjvVbK2TC6mrYV+dHUG/FZoMWVqWrHa4F2bradrPRA7pWm\n/j6nBtrYNWFq5h7Vk08tm1kC+2PwKK3jOGWYTLzCc4fnOI5JJm5L8a5ljuOY7FWCx46ISF0ReV9E\nZovI5yJy7Q7P9xGRrSJSrZB2a9QUbK6ItC6kNxWRWdFz9jpJMXGH5ziOyZ4leBhsAq5T1SNJFHW7\nSkSOgIQzBM4A/ptvLCINSdTVbAi0AR4Xkfw1wcFAD1VtQKIOZ5tdnVNagxZPa5haZC1YP04v8xgP\nSbiAPET7mLbZvG/qw7g40JaZOUXwT/5s6pfyTKA1wg44nEnYi2AWjUzbPq8/buojzgnr03XcYPco\nuLXS3aZew6jB9zEnGpZQj68CzQq+APQRu3ZhdpidxqAn7QYRf97yhKn/6tkw9Wr2ZXZKl/UZtl76\ntmmbWztM1QNo+dTkQJNG9v+Hbs0Hm/qQleF3t1GN8Lhg9/JYOSgm6BRmYCaoFEpPduoSaFfI8KSD\nFqtLsOBVfXPRQQsReQV4RFXfFZEXgbuAV4GmqrpGRG4FtqrqwMj+TSCHhFN8T1XzneUfgWxVtf+z\n7gRfw3Mcx2TPkniHzfFPRZXRmwCTRKQ9sERVPyu4gAMSmc2fFPp9CYlmYZvYvn3E0kjfJdzhOY5j\nsldW/HMfbkk8doaI7Ae8BPQGtgJ9SdzObjNJ5hxLijs8x3FMirrCa7UntCr0+91rQxsR2Qt4GRiu\nqq+ISCOgHjAzurqrA0wTkeaEjcHqkLiyWxr9XFhfWuLJRLjDcxzHZK+YPjHFIQo4PA3kqepDAKo6\nCwoWMkXkKwrW8MYCI0TkARK3rA2AyaqqIvJ95BQnA10hpnZZMUjK4YnIIuB7YAuwSVXDjtWO45RN\nkrscOgnoAnwmItMjrW/UCyefbREiVc0TkVFAHokVwV5aEFHtRaJSYUVgnKq+uasnlVSUtrCHNp7T\nuzXsfGV1Wm+17D/m8WfXCiN038Z03rqb20z9jZnnBtqFxww1bS/hWVOfahS9vHaL/UdmbFbY+eyi\ntXYn+mVV7KKl+xppTHEpT/dyg6n/mTASWp3Vpm1NVgRaXBT68Lv/a+pZPdYHWueaI0zbCmw09afe\nuybQ+rYKC5kCzKBxoI3rf55p2/cO+xjjODvQVhtd5wDmbjzc1G+scG+g1Ypp13b7pPsDLate+L4B\nbHlvP1Pv3imM7Fvd6F6VzklHabUE7Svlf+UntewXP0nHcXaBDFzwSnbjsQLviMhUEbEz8h3HKZtk\nleBRRkjWh5+kql+LyG+At0VkrqpuK1vyTk7Breoh2XU4pLRbyjtOOWFVbh6rcvNSe9AMvMJLakqq\n+nX07zciMgY4DgrqNJ2ec0JyZ+c4TrHYP7sh+2cXlKWa12908gdNIkr7S2WXb2lFZF8RqRz9XAlo\nDTEr3Y7jlD2STKb9JbLLUVoRORgYE/26J/C8qg4o9Ly+oWEO4y3cE2gfbTzJHOPPFcIimw2xL9v7\nTnvQ1Ls0fTLQro3ZxvMWZ5q61TLvfF4ybVsPDAuRHnezXcRy8kw7x3PgMWG0sqaRGwvwQ0wRnyuf\nCCPRU3raLR2t+Z270i5MeWaNsFgo2Lmtce9RNrmmPolwV9PLnG/aTv5vmBd81kFhgU2Ajtu+pttj\ntTw8a559bjF1Zxl2Tnh+cVFaK/LdYOMC0/a6Cvb3+VDmBdqP7BtoPWRk8lHa40tg/0mGR2lV9Ssw\n9gY4jpMZlKFgRHEpQxejjuOUKhnoHTJwSo7jpIQM9A4ZOCXHcVJCBnqHtE6pzYJwsf7BBmEa037d\n7TozV458oNhjzW9q7/F7gp6B9vu5n5u23Q63U87yFjcNxQ32eay+OQwiLOR3pu1Lx9gL8lbwpB93\nmrYtPvjU1L+/LCy8HVeItOvakYHWoYadFhbX8W09lQOtP3eYtjdsvM/Ub6/wt0AbTZgaCDD9oHD5\n+AX+aNpmYX+/LtjwYigeHhZDBRih15n60SXYmPCbx8I0souuesG0vbOr3WJu+bAqgWZ1AkwJGbgt\nJQN9uOM4KSEDvUMGTslxnJTgUVrHccoNGegdMnBKjuOkhAz0Dhk4JcdxUoLf0paM2Q3CAp7rjGje\n+iH2OzvdSORo8ZUdldyu31Ehnu/UOdDmHW5HtfLGGNFY4OyOYQHPGkbRTIDHNoatCavuY0eFT1e7\nmOnAxTmBVrfufNO2bUs7ncoqbvkJdq7QkCphm78T+di0PejOb0ydk0K9X2s7snxfBbto6YTFrQPt\nybrhuQHszc+BNuzFK0zbGhfYRUvXrakRinNMU7KxU+1yPzorFGO+i8P7hAVKb8COWHcYZkfJX7k5\n/D7fN7CtYWm3iiwRGXg55I24Hcex2acEjx0QkWdEZIWIzNpBv0ZE5ojI5yIysJB+q4gsEJG5ItK6\nkN5URGZFzz2c7JTc4TmOY5NcAdBngTaFBRE5FWgHHK2qR0Hi8lZEGgIXAQ2j1zwuBU1rBwM9VLUB\n0EBEtjtmSXGH5ziOTRLloaJCwN/uIF8JDFDVTZFN/jpIe2Ckqm5S1UUkatM0F5EDgcqqmn9//hzQ\nIZkpucNzHMcm9fXwGgCniMgnIpIrIvndsWqR6EGbzxISrRp31JdG+i6T1mXJSTQPtK+oF2j7zbBT\nf65oHt6yt7jQDlr0m3KTqT9GGER4EDtNiC9t+Y0XwvSm1y46zbSt9MHWQGtgZ1hx6Fth7TwAbg+l\n3CnhewlwN31NfcXacEH+pCofmbZWvbi4xfR+/ez3+botYf22uPe5M/aCPEZpuGV1wzp7AG0xgjUx\nqVDfzI5pv2XEgfp1tOf3ZExKHUapw0f79jBNr1rzdKD9UMm+5qhcwe5mhlE68v77jS8Mdqe2ElFE\nlDb3v4lHCdkT+LWqHi8ivwdGAWFkM41kYBzGcZyUUIR3yK6feOTTb2KxjrgEGA2gqlNEZKuI7E/i\nyq1wMnydyHZp9HNhfWmxRorBb2kdx7FJ/S3tK0ArABE5FNhbVVcBY4E/isjeUSX1BsBkVV0OfC8i\nzaMgRtfoGElNyXEcJySJaikiMhJoCVQXkcXAHcAzwDPRVpWfgYsBVDVPREYBecBmoJcW9J7oBQwB\nKgLjVPXNXT8rd3iO48SRhHdQ1U4xT3WNsb8buNvQp0FMXbNdYKe3tNYGQhGpJiJvi8h8ERkvInbK\ngOM4ZZfy2LVMRE4G1gPPqWqjSBsErFLVQSJyM4nIyy07vE6XaLXgeLWXrQm0xrX+E2gAn2wI+9ru\nc799nhfdMcTUv6B+oE075GTT9tEv7eja1f3D6NrMOxqYto1nh6G/K460N4ivwEhtAl59N/zj2O20\nwabttTxi6sOMP6QPfm5HdP961K2BNj6mg9vImCKbh/xleaBtGGj/PW1SwY60L3j+mEDTwXYjrDcn\nhh3f2lxud4eLi74//G6YiraCmqbtjzHd4brzbKCdG9Ml7WQ+DLShf73SPrmptrz1+fD9+LZamOpQ\nXX5KvmvZoyWwv7psdC3b6RVezAbCdkB+eeChJLkZ0HGcXyAZeIW3q6daU1Xzs+dXQMyfRcdxyi5l\nyJEVl6SnpKoqIuZ98f05P2z7+YTsvTgxO+yz4DhO8kzM3cJHueGm96Tw8lDbWCEiB6jq8ijfbaVl\n1Ccn7IjuOE7qaZGdRYvsAg81qJ+dvVQijCooZZ1d3Xg8FugW/dyNJDcDOo7zC6Q8ruEV2kC4f6EN\nhPcAo0SkB7AIuNB67QtcFGgn1goLS1blO3PsRZUOCrTn7rjYtB31RDdTb9BzZqB98+V+pm1c1LRi\n79WB9uMT1U3b9j3DlofvcLppO3vDUaY+77SwcOYrdDRtj51rV6xsOjSMvLYdYBcLvcSINC7iYNO2\n5kbzYp6HHwojno2ZYdoumBdGYwGwNjfFRAqPY1Ioxm2UvcyWez/2f4EmP8XsWrA/Qv51TPgd/7pm\nuDMAACP9qtZdX5umE7B3EvyD4r7PKSgAWh5vaYvYQBjzFXAcJyMoQ1duxSUDp+Q4TkrIQO+QgVNy\nHCcllMdbWsdxyikZGKVNq8O7fnyYDvVk67AL1c/sbb7eSvOZz6Gm7aieVucmeMtIkToWO7XpDvqb\n+k9zwxQ5Odxe3J5jFDh9KmbV/Df2bh4zfevOFweZtrMvsOsnrh4QpkLdx42mbTMjj6kCG03bpyrY\nc7E6zNVimWmrq+wMpKXnhO9zh5hioXVZHGiPnmEX7/xL+wGm/uCdYaqdfmmf29A+ZlyObmNGBdrg\nFXYArdfCMDg0H7s46dm8burDCIN293CzYRkTOCkJfoXnOE65IQO9QwZOyXGclJCB3iEDp+Q4TkrI\nQO+QgVNyHCclZOAanve0cBzHJsnUMhG5VURmi8gsERkhIhWKKh4c2S8Qkbki0jodU9ppAdBdPrCI\nwieBPk0vCbRZMRWcL1/9ZKBtGv4r01Ybx0TXWobRtfu4wbTtyROmftmGsABot0pDTNsXunYPNLnF\nfo+HHXm+qXcd+lKgTejW1LQ9eeg0U+enUFrU006ds97/IYSfE8DQjXYEcr/hRrJ6tj3vW+uHqXMA\nl8hdgZalB5i2hzQJC45+Od22rbs2tAVYVKVOoDX4YIlhCaNbnmXqbau+EWh7zTIMgaZ1w7aczWNS\nwF7YEqasAay5x2jL+jfD8CdJvgBoTOFU0/6Q7QuAikg94D3gCFXdKCIvAOOAIzGKB4tIQ2AE8HsS\nvWffAQ5V1ZSWgPErPMdxbJK7wvse2ATsKyJ7AvsCy4gvHtweGKmqm1R1EfAFcFyqp+QOz3EcmyQc\nnqquAe4H/kfC0X2nqm8TXzy4FoletPksIXGll1Lc4TmOY5OEwxOR+sBfgHoknNl+IrJd1kHUirGo\nNbWUr7d5lNZxHBMtIkqbOxFyPyry5c2Aj1V1NYCIjAZOAJbHFA9eCtQt9Po6kZZS0hq0WLs5LOl+\nVdZjgTZsXFjjC+CnsDEV51Z62bTtGNMp6vLhwwNtdBd7AfrcGeECNMAnjcP6bXGBDysdbsLv7YBT\n0ynhIjbAtJVhLbTuNR43bdth17g7d2o4F/nA/qzn96kbaHvHpJb9ds03pr7Z+M9xZBW7Ht6JhDUR\nAebQMNAmjc82bakSSrJnzHe5akwaYP2w5l+9Df81bSu+Yh/j6z+FRfxqrllr2uZVC9MA11HZtD1h\noP3ezbg5TK28iXsDbbx0SDposcmehsleVYKgxTHA8ySCED+RaKY9GTgIWK2qA0XkFqDqDkGL4ygI\nWvxOU+yg/ArPcRyTLck14p4pIs+RaDi5FfgU+D+gMkbxYFXNE5FRQB6wGeiVamcH7vAcx4lhYwW7\nqIfNz4GiqoOAHaterCGmeLCq3g3cXYJBS4w7PMdxTLZkZV6qhTs8x3FMtmRgbtlOt6WIyDMiskJE\nZhXSckRkiYhMjx5t0nuajuOUNpvJKvajrFCcK7xngUeA5wppCjygqg8U9cILsl4MtLf6dwi0Q++w\nI1LPSFhUcvT680zblyrZesW2Ycexa3nEtP2usdU2Cx7kukD7aMtJpq1F1Yl2Z6onudzUG9f4T6Bl\nsdm0HbKgl6n3axYWw1zUzE4tszqiPUdX07ZRNTtvasjn4Xn8qYpdvHMZtUx9EEYBz5gVnea5uYH2\nMmebtufWtqPv1ywdGGiLKtnd2o79k9FyDLiUZwKtZzU7RbHDzLdC8SnTlAmP2KmE7RkbaAO4JdDG\n24ctEVsy8AawOF3LJkR5cTuyyyFvx3F++ZTLW9oiuEZEZorI04UrHjiOkxlsIavYj7LCrl6zDoZt\nDSDuIpEz12NHoy9yCppSV8s+imrZdlUUx3GSIy93FXm5q1J6zI0xvWbKMrvk8FR1W/cZEXkK7O3+\nv8uJ6+HtOE4qaZi9Pw2z99/2++h+85I+Zrlcw7MQkQNVNX8lviMQUwHMcZyySlm6VS0uO82lFZGR\nQEtgfxLlXO4EsoHGJKK1XwE9C5V8yX+dzjSOd8yQcLy41nhsCKW/3dfHNI3Lz2xIXqD15mHT9kyM\nKBpwMmHOqxW5BXj87vD8Xuz7B9P2wrZ2Huxrr50WaFZbQohvAfnIpLB1n+xjf9YTjgkjgj+yr2nb\nKOZv2wTC/N8bjRxPgNGca+rP0j3Qjo4Z7ywJc6RbqR3tX7jid6a+ZcB+gXb+Q8NM2zh6EeY4Hyph\nlB3gfg1Dzh14xbQ9Y/Xbpt6wevh9fpELjHNYknQu7ad6RLHtj5U5SY1XWhQnSmvdl4axeMdxMoqy\ntL+uuGTeTbrjOCnB1/Acxyk3ZOIanjs8x3FMfvZtKSXjaCO+MKxb2KlrEFebr7cW5OetCdPNAPb4\n2l6QP+rIKYFWnTDdDODSziNN/ewRYdHRcbXtVLbK81cG2rpldkqXdo5Z4/0glP7dspVpehjz7WMY\nU5n/sD3eQg0rrZ621jgJ4Nwqo0w9zyje2Rk7tazpB+HCO8D648MrikpX202rZFFYjug17JTuPzzx\nnqljFJiNC8r8QXascpTgYN0n0O4yghMAbY3dW3nYgYEK+9gFWGfUPiHQ6i61vgNhodCS4mt4juOU\nG3wNz3GcckMmruF51zLHcUySzaUVkTYiMldEFkRNt3c7foXnOI5JMmt4IpIFPEqinPtSYIqIjFXV\nOSk6vV3CHZ7jOCY/UyGZlx8HfKGqiwBE5F9AeyBzHV79+z4PtIu3qyOaIC6l66YXHw20By640rTd\nGlOe73nCaGqXQ+1Wj6/Ot9spjnsvPMbspWHLPYB1/cOI7L/vsCOsL3UKI9YAJxlpco2x06YGGsUf\nAa65IKwsecBDX5q2y68J53LlI3Zt13u50dQPnxS2NxzU3I6+X9hyqKlXmhpGZC96cohpqw+HWyay\ne9uFPsfdcY6pWyyinqnfGdYKTTD8p0A6r8tLpmmLFz8NtNMOtVMidbr9ff5g6XGBdplRRTTHfHXJ\nSHINrzZslw+5BGie1AmlAL/CcxzHJMltKelpeJ0k7vAcxzEpalvK3NwVzMtdEfs8iXW7wh3e65K4\nytutuMNzHMekqFvaBtm1aJBd0JtkbL9g+Woq0CBqD7EMuAjY7QUy3eE5jmOSzBqeqm4WkauBt4As\n4OndHaGFNDs8qyvX3oQpM83XWJXzQMaHywD7tF5jD1bNlq2Upy7v2kGLhdS3z6O2sRwRs4g95456\ngVaTMN0M4LGN9qL+wRW+CrSlaw8ybW+ocp+pT6gQBmCkrb2s8tBrPQOtd+f/M20fHnGFqX/cfF2g\n3fy03R1uUQ871e7hZuGxrXqGAHJ8OJcfN9gL/e9UsoNGVre2cTfYKYNWAA7sbnLz+9vpjxhxrpuO\n6Weayr/sz+pYwu5pb3NGoOXYZ1Aikt14rKpvAHYkaTfhV3iO45hsTG5byi8Sd3iO45hkYmqZOzzH\ncUzc4TmOU27IxPJQRRYPEJG6IvK+iMwWkc9F5NpIryYib4vIfBEZ7424HSfz2MKexX6UFXZ2ppuA\n61R1hojsB0wTkbeBS4C3VXVQVAXhluixHf25IzjgvvwQaG8MtrtYaW8j6mY3HKPhHdNMPW9G2JEL\nO0jIs1xi6qsPqxhow26+2LQ9fFyYYiWH2UUs69T/wtRP553wuFXsKGE7xpp6+zXjA+2lf9tRzHWW\n/K5pyl9efcLUH2kfFmvVtvZ4y6li6lZaVxZb7BPZ34jg228Ff6hpFwC9vNWTgaZ9rzFtH+B6U2/f\nInyfGW2fx4Aafwm06tjNs/XHmOKwxldUOlgR3eQbiJW7W1pVXQ4sj35eLyJzSOTItaOgXuxQIBfD\n4TmOU3Ypdw6vMNGO6SbAJKBmoT60K4CaKT8zx3F2KxvLa0+L6Hb2ZaC3qq4TKbhcVlUVEXOX5KKc\ngmbJVbOPpmr20cmdreM4NrNy4fPclB6yLK3NFZedzkhE9iLh7Iapan6b9BUicoCqLheRA8FOJaiX\n0yV1Z+o4TjyNshOPfF6wMzhKQibe0u4sSivA00Ceqj5U6KmxQLfo527AKzu+1nGcsk2yJd5/iezs\nCu8koAvwmYhMj7RbgXuAUSLSA1gEXGi92Cr2eTbjAu3AN+0orZwf3ilrMzv6lDfGiMYCZ3cM82bj\n2gfO+iAsrggwsmX7QOu92M41/fLsAwLtSh40bQe/bkf+RjwettiTdnZuZV5PO9I7oHUYEZyk9vxe\nlT8G2iWz7PdZj7B1eTU8v2XtaxmW8dRnYaBdPeZp0/bMjkbRWDvgSdPmE0z93xLmTl+udv7voJjC\np/UnhtHOuSIxAAAdI0lEQVTzF7nAtJ1ntE4cck4v01bGhf9PAOro7wJNh4efSfIx2szch7ezKO1E\n4q8CT0/96TiO80uhXK7hOY5TPilLt6rFxR2e4zgmP5fXbSmO45Q/MnENT1TT02tDRPQR7RHodbdr\nZJTgY040jzFIjgq0afpX07YyYQFKgAZfGWX0jWwggFd72l3LrPOrynembS7ZgfbiFnsR+/mszqZu\ndSKrjx2cqMDPpm4FZrpcYhc+lSF3Btps7G0ND28XrC+gEbMCrRlTTdvjR9oFX5kbSg/0s7vUXf/U\n4EAbfdlZpm3cBtre/CPQptPEtK09JqbwbPgV5c0GLUMR+Mj4Hh3MItP28hVh2hvAzTXDyrPWrecg\nyUFVdzl2ISJ6ng7fuWHEy9IlqfFKC7/CcxzHJBPX8Irch+c4TvklXfvwROQ4EZksItNFZIqI/L7Q\nc7eKyAIRmSsirQvpTUVkVvRcTAmRneMOz3Eck81kFftRQgYBf1XVJsAd0e+ISEMS3c0aAm2Ax6Ug\nj3Uw0ENVG5DohtZmV+bkDs9xHJM01sP7GrbVCKtKooctQHtgpKpuUtVFwBdA8yh9tbKqTo7sngM6\n7MqcfA3PcRyTNG5LuQWYKCL3kbjoOiHSawGfFLJbQqIc3Sa2b+K9NNJLTFod3iSaB9o1TZ4KtGOn\nh63nAHT+yYEme9lR5R+/iwkQLQil7J5257hLeNbUBz0RRjHv6nmDaftm9fAPzx7X2+c84TY7HW7R\nOUeEohENBODMGN3IvMp+1p63HhVGN7/ps59p+8SQMGUNYGb3BoF2NY+ZthP3CNsKAvB9KFXbsjQU\ngevXhFHapwiLkAKM+8Buvdh5SZgCfs+fbjZtH612k6lri/B79+0KuwD43ybdHWgPNQ9bZAI0rznJ\n1K0iqYcxz7RNlqJuVb/Lncna3JhoOxAVCg7zLOE24FrgWlUdIyIXAM+A0WsyDfgVnuM4JkXdqlbO\nbkrl7II/2Iv7bb+FRVVjHZiIDFfV/NTUl4D8q6ClQN1CpnVIXNktjX4urNt/CXeCr+E5jmOSxmop\nX4hI/mbFVsD86OexwB9FZG8RORhoAEyOKq9/LyLNoyBGV3axQpNf4TmOY5LGfXhXAI+JSAXgx+h3\nVDVPREYBecBmoJcWZEb0AoYAFYFxqvrmrgzsDs9xHJN0OTxVnQrGAn/iubuBYLFTVacBjZIdO62p\nZXdruMD9LN0DbeqW3wcawK/Gbgq0iR2PNW2/xq69treRerV4u2WCAl6jramPfzWsh7e6fdjJDOw0\nNCvtCuA77MXtY8aGkZYO7ewafmPW2OlpTar9J9BmPna8aXvIVbMD7Ys1dpTEOi7YaWRP1bY7gLVd\nOsrUO/N8oHU691XT9oPRYW2/uKji89jv0bWEte+OfW+OaSv17RS+uw66NdCsrnMAeTQMtBP52LS9\nMya1r6LR9e9G7gu0o+TLpFPLjtBPi20/R4711DLHccoumZha5g7PcRwTd3iO45QbMrE81M6a+NQV\nkfdFZLaIfC4i10Z6jogsiZJ/p+9qXpvjOL9c0phattvY2ZluAq5T1RlRb9pp0Q5qBR5Q1QfSfoaO\n4+wWyt0tbbThb3n083oRmUNBDttOIzJWsc+PaBFov5oTRmMB5PswgnwEduQo7vL7FsKCiS9xvmn7\nwUf2hapODqfatL3dCasy6wPtn/zZtK1pt/OFMJONJu1mmKZ51Q4x9QpsDLTfXPU/03bh5WFEduaT\nYaoYwKm8b+rHEaZCyU32DoA+/M3UX+CiQOu8LmZ/6buGdvhPtm32PqZ80oIwQprXKoykAugQOwI8\npXv43j3IdabtyKVhtFjvtaP9o9p2M/VWp/070MbSzrC0C7WWhHLn8AojIvWAJiSSe08CrhGRi4Gp\nQB9VtUsAO45TJtn4c+b1tChWall0O/sS0FtV15OoTXUw0JhEqZf703aGjuPsFrZs3rPYj7LCTs9U\nRPYCXgaGq+orAKq6stDzTwGvWa8dnZO37ecjsn/DEdm/SfZ8Hccx+DJ3MV/mGv1bkmDL5nJ2Sxsl\n6j4N5KkWdG8RkQNV9evo145gpxKcm2OvhziOk1oOya7LIdkFGUTv9bPLS5WEcufwSKzVdQE+E5Hp\nkdYX6CQijUlEa78C7KJejuOUWTZvyjyHl9ZcWh1mPGGkaNZtPD8UgTcJo6bfxuSf/jqmbeIP7Bto\ncS0dv6C+qXfa8K9Ae7uSXe7r+IfDooif9D7Gtr0wpoDiSaF0cm+7t+TEjvZ5DBkTRjx/wI4IDuPi\nQFtFddO2GdNMfcTdYUvOwX3tSOOd9Df1f/HHQItrTznHyEs9gjzDEup1taPhPz8aRt9fqfIH0/aC\nmmF0FLBXryvZplZ92b5j7bajd194l6lXfjacy7ohNQJNribpXFqWxkS9LWrv47m0juOUYcrhLa3j\nOOWVnzLPPWTejBzHSQ2bd/cJpB53eI7j2LjDKxmDulwdaA2NhWWrSGecPpBbTNvXWlxon8SjoSQz\n7UBNl25PmvqBlZYFmlWIEaBp7zDlbNr4sPsaQL9Rdies/ivuCLTvN/7KtK3UYqupz+WgQKu34b+m\n7bpKlQPtIl4wbWtu+MbUMWqLXrJhqGl65Yu2Pqh7+H35O31N26ONnVBt3vvAPje70RrTqoTBpHuw\nu5a9v+JUU7c63cWlLg7KC3MGx2N/NwaNCt8LgA5GKwdpYX2fUxA/SJPDizqV5QCHA79XLag0KiJH\nA08AlYGtQDNV/VlEmpIo8b4PiRLvvSP7CiT61B4LrAYuUlX7i4438XEcJ45NJXiUjFkk9u9+WFgU\nkT2BYcAVqnoU0JICtzsY6KGqDYAGhSo09QBWR/qDYCTPF8IdnuM4NltK8CgBqjpXVa29aK2Bz1R1\nVmT3rapuFZEDgcqqOjmyew7IbwDdDsi/ZXgZOK2osd3hOY5js7kEj9TQAFAReVNEponIjZFem0R/\n2nyWUlC1qTYkyjKp6mZgrYhUixvAgxaO49gUte94Zi58lhv7dFQ38wDjqb6qaubeA3sBLYBmJNo3\nvisi04C1xTnd4uAOz3Ecm6Ku3I7MTjzyGb59lzVVtVOAimYx8KGqrgEQkXEkghHDgTqF7OpQcMW3\nFPgtsCxaA6yS/3qLtDq8yUbryR/ECJvGZLc1uDCs/tD0RTsa21djUnT6hyk6eq4dwdqw0b7Dn1qh\nWaDdRRhJBfhwY8tA6976cdP2sY12JK5nzScCrRt2ZHNZH7s9ZV/+HmjVK602bW8aFH4mcqT9obx8\nztmmvrFVWDttHOeYttd1f9DUrVLhz9PFtLVaHv5kdjoFXWl/3p9yRKD9joWmbVxa3nFDjboZ3WOq\nlhito+OKhdZghanXItwxoNXC+aUkx6t0tqUUPtW3gJtEpCKJUEhLEpXVl4vI9yLSHJgMdAX+Eb1m\nLNCNRJ3O87FLw27Dr/Acx7FJ37aUjiQc1v7A6yIyXVXPUtXvROQBYAqJy6DXVfWN6GW9SGxLqUhi\nW0r+n4+ngWEisoDEtpQwIbsQ7vAcx7Ep+XaTYqGqY4AxMc89D2FHdlWdBjQy9I1AzCbcEHd4juPY\nlHC7SVnAHZ7jODYZmFqW1np4b2i4gN9mZZj+06uG3RLjTukTaAdMjTnfsGQdACPu7RBo93GDafsz\nFUz9dN4JtIdev9W01f2M5WK74Rj1e39u6h8ZBfEOnGfX+9Mp9vL03meFkfxG1T8zbd8mDKhV2Gin\n++33L/vPfp9uYScyK7AAMGmLHV2omhXO0XovAA58N7R947Rs07bNi3bK2Q0XhAGtqjF1FbNiLne+\nM+ozDlpoB7T0szAo1rOj3V2sudEFDuDS/iND0WiKJzVTUA/viRL4hp7i9fAcxynDZOAVnjs8x3Fs\n3OE5jlNuyECHV2QurYjsIyKTRGSGiOSJyIBIryYib4vIfBEZLyJ2ownHccou6auWstso0uGp6k/A\nqaraGDgaOFVEWgC3AG+r6qEkdjbbReocxym7pKlayu5kp7e0qppf6XJvIAv4lkRJlvwQ7FAgF8Pp\nnfV6bnC8Uee0DbQmTA80gNl6QqANiymu2GX9y6bellGBNq2RXXTxmVmdTP0r6oWilRYNfNN0v0B7\ntmV30/ZDTjH1x+gVaFce9oBpe9phdh72z4dWCbRP5tvd01rxfqDNGB2+9wD805bvax6m9jU83O5w\n9nzWn0zdKvi6DDt1jnphBPF97CKdbTbaUdoJRvHN7lZrMeDK2+zUvvP/Hrbm0xz7OqLnsDAi+8Sg\nmOqkZqdneHNYuPPhdax0P7uQaYkoQdOyssJOy0OJyB4iMgNYAbyvqrOBmqqan+y3AqiZxnN0HGd3\nUPrlodJOca7wtgKNRaQK8JaInLrD8yoi6dnM5zjO7qMMrc0Vl2JHaVV1rYi8DjQFVojIAVEVgwMB\nu9Px8zkFPzfKhqOzkzhVx3HiWJK7kKW5X6b2oGVoba64FOnwRGR/YHNUxaAicAbQj4KSLAOjf8PO\nIgB/yknluTqOE0Od7PrUya6/7fcp/cLsoBJThm5Vi8vOrvAOBIaKyB4k1vuGqeq7IjIdGCUiPYBF\nlKBageM4ZYTy5vCiZhrHGvoa4PSdHXyPZhsC7cKHw6hitauXmq+/OOu5QHtoob0DpstKO/7S0ahC\nIzfYS446NSYV0MiFffcye/o1Xl8XaLefY7ca/Cc9Tf2uhwcEmtX+EeBPjDD1DbPC9+P4F2eatj9c\nEBa3rG4UXwVo2cmoYgnQP5TmVAq+OgBM7RMWVAW4fVCYUy1v2J/V1++HWz/XEbabBDuvGGDVlv0D\nbXNWlmkrK+3zuN+IqHO3aUp1wgKsw286z7QdQWdTX031QJt0Z3agGWV2S055XsNzHKecsXF3n0Dq\ncYfnOI5NBt7SeptGx3Fs0pRaJiL3isgcEZkpIqOjLW+Fn/+tiKwXKagPJyJNRWSWiCwQkYcL6RVE\n5IVI/0REDipqbHd4juPYpC+1bDxwpKoeA8wHdiwu+QDw+g7aYKCHqjYAGohIm0jvAayO9AdJ7ByJ\nJa0FQGdog0CvYKQPzQpL1QNwnFEEcUJMOla7LWNNvUqHcLwLX7PThFYZC8IA7936h0B7eoC9qHwy\nYXChMmEgA+ILkVrnEVcQ0no/AabTONCu5RHTdjF1A20FNUzbgTFp082YGmh2yhMs73eIqe919feB\n9mT1y03bd4yY2cWEQS6AD40UMrC7pN09JiwKClD17K9N/btZBwbaRc2GmLb1WBRoT225zLRdkPU7\nUx9BmJb3GmG65njpkHwB0LYl8A2v7VoB0Kihz3mq2iX6vQNwIrABWK+q90d7fd9T1SMimz8C2ar6\nZxF5E7hTVSdFbRq/VtXfxI3nV3iO49iUTmrZpcA4ABHZD7gJyNnBpjYFfWgh0Yu2dqHnFgOo6mZg\nrYhUixvMgxaO49gksS1FRN7GLrHRV1Vfi2xuA35W1fy9VTnAg6r6g4ikpVy8OzzHcWyK2pbyXS6s\nzY19WlXDRimFEJHuwNnAaYXk44DzRGQQUBXYKiI/AqOBOoXs6lBwxbcU+C2wLLqlrRLtEzZxh+c4\njk1Rt6r7ZSce+fyvX7EPGwUcbgRaRjU3AVDVUwrZ3AmsU9XHo9+/F5HmwGSgK4lG3lCQ5voJcD6J\n+pyxuMNzHMcmfZkWj5Cor/l2dOf6H1U1Ula2oxcwBKgIjFPV/JSfp4FhIrIAWA38saiDpDVKu9VY\nOvzH6isCbQt2Os9Go21i34EP2gNm27L+z1gKGGfb9n02LGIJ0IiwveErdDRtD2VeoMVFR+PSh/bl\nh0C7fMZw07bqEXb08JufwujhtCp2AVCrNeGzXGLaDmySY+ojp7cPtLc407Stm1hjDjjd+ON87pbR\npu1hWeH73Ba7GGpWzKWKlZZXq7p9N3TPartQp1W81tpdADDEeE9XxJSSPJkPTX0hYfS295j/CzQ5\nNwVtGpuUwDdM9zaNjuOUZTIw08IdnuM4Nu7wHMcpN3i1FMdxyg1eLaVkdFwdLgrPoEmgfbblaPP1\nN2bdG2hn3WwvYtfnC1M/u3nYzeyNleeatifwnqmfyVuBtsjqZAZcRxhUOYmJpu3KjfaC9ZsV2gSa\nNNlq2hJTw2+fOusDbXgVu1tYp89fDbS9j4r5ttuN3ahpVPmPSxk8xUi/A/s9zcqyEzU7GwGHq5c9\nbdr2qhXW2QOoPdcIUIQfNQDTaGrqPY02bqsJ6+wBXD91cKANbtbNtG27dMdU0gRX1H4i0Lp2DIMW\nEAYHS4zf0jqOU27wW1rHccoN5a2Jj+M45ZgMvKUtslqKiOwjIpNEZIaI5InIgEjPEZElIjI9eoSL\nTo7jlG3KWyNuVf1JRE6NqhfsCUwUkRaAAg+o6gOlcpaO45Q+5XENT1Xz85z2BrKAb6Pfd5pG8krn\nMHVq0IirA63m+uXm638cFBbCvObvdkHTSxhi6o8bXaW0ln3qMs1OpVnZNCyGOemxbPsYj4bHGDGn\ng2n7Y4V9Tf34MWF3sTf0VNPWSr8DWEj9QLtyYxhRBOhc+cVA6xNTOFam2+/R1QwKtGmL7cKbdevO\nN/VaLAu0q3jMHm9kGJE9tJPRXg64jKdM/dXDWwfahJhiodON3QUAPxJ+hg1m2B3fMC4PTh9h57rr\n7WEnOQDpYrz/dULJo7Q2Oy0AKiJ7iMgMYAXwvqrOjp66JqpJ/7SIhD3zHMdxfmHs1OGp6lZVbUzi\n78gpIpJNor78wUBj4GvA3ujkOI7zC6LYUVpVXSsirwPNVDU3XxeRp8AuU5Ezq+Dn7BqQbe+zdRwn\nWSblwuTc3X0Wv3iKdHgisj+wWVW/E5GKwBlAPxE5QFXzF946ArOs1+fYG+0dx0k1zbMTj3weK35B\nzngyL2qxsyu8A4GhIrIHidvfYar6rog8JyKNSURrvwJ6pvk8HccpdTIvapHWAqA3aU6g/87Iea3P\nQvMYrQ75Tyj2t8cb0MUu0GhF3fryd9P2O+zYyx/eC3NsKzcPc0cB/lopPEGrDSLAI8tuNnVGhtLy\nPlVCEWiwYYGpn1kpTAq18k8Bzns1rIg6rX1D0/ZyrLxNu8hp3LxvW2Mv+e4xJ/wuzjjpUNP2MGPe\nP1Tax7S9lGdM/WKGBVrcOff+sz3vy/4ZzjsuKnx8/zD63v2Ox03bIYvtAsAL6oYh2euN8O+/5cLk\nC4CytgSvqOIFQB3HKcv8uLtPIOV4X1rHcWLYVIJH8RGRu6ItbTNE5F0RqRvpZ4jIVBH5LPr31EKv\naSois0RkgYg8XEivICIvRPonInJQUWO7w3McJ4a05ZYNUtVjou1urwB3Rvo3wB9U9WgSncgKrzkM\nBnqoagOgQaF01h7A6kh/EGJ2zEe4w3McJ4b0XOGp6rpCv+4HrIr0GYV2f+QBFUVkLxE5EKisqpOj\n554D8tOX2gFDo59fZvs+twHpLQDKmEA7fni4cNu4ixGcAOp9+VWgLeJg0/Zn9jb1oYQFFuO6afVf\nfYepr2oZFnR8IOt60/Z5wiKb93KjaXt2rbA4KQB9QmncE+eZpr162oveAyflBNrg5naxST0oXGue\nwlGm7acntzD1bh2MlKwwMxCAet3nmPrkk8J9TFaXLoDrK4UL9XHphc9tsed9fVZ4jLhUvU0D7aDF\nU4uvCbQNNezriJPvGB9oE8+y+1Xv+YZdm8kKqjQmTKn7t/nqkpK+KK2I/J1Ef9kfgOMNk/OAaaq6\nSURqU9B4GxLNt2tHP9eGRBs8Vd0sImtFpFpcM24PWjiOE0NRV25TgKmxz4rI28ABxlN9VfU1Vb0N\nuE1EbiFxK3pJodceCdxDYt9vSnGH5zhODEVFaY+KHvlsX5hCVYvrrEZQqFO0iNQBRgNdVTX/Fm8p\n25dIqEPBFd9S4LfAsqiiU5W4qzvwNTzHcWJJT9BCRBoU+rU9JLqZR0VIXgduVtVt61yq+jXwvYg0\nFxEhcSuc34hlLGxbtzofjG7uhfArPMdxYkhbatkAETmMRBH5hcCVkX41UB+4U0TyI7dnqOoqoBcw\nBKgIjFPVN6PnnwaGicgCYDXwx6IGdofnOE4M6QlaqOr5MfrfgL/FPDcNwjZ4qroRuLC4Y6c1texz\nPSTQx9Ax0G7/yk41Ov/gMPVn8LY/BttzHzeYutU28cBff2fajvjWLtTZ6cWwjaH8yn7fdFoY8ZQK\ntu0jfS4z9alGS8C41obPExZZBfg7twdaVex5H/HuokC78jS7mHVc+t0PhAUrb4nZEnXCNLtQp44N\n37uH+9mFLI+RMGraV+27mf9c08oe74xwvDPbvWLavjXG/m5kd3wj0KzvHEA9wl0HZ/COaWvtcAB4\nbm3XQHuxygWB1lbeS0FqWZhyGM/ZnlrmOE5ZJvOKB7jDcxwnhvJXHspxnHJL5hUPcIfnOE4MmXeF\nl9agBVPDY2uWsa55uX2MiVOODbQWX31q2jY/ONfUOxAuQvddOsC01XV2p6hBh4ed1q7aaKd03Vkh\nrDS7jsqm7RNd7Rp+1wwLF/sfWWnXzmtYY5qpL9tYK9C6VnjOtH1krHFsuwwdMw9vYOqPcG2gxQVJ\narDC1LcYf39vXfCQaftMg06BdukCo5AgUKPBf039xw3h5718Tys5ACrN2GrqU5qHKXhjaWva9l8T\nfu/mxKTfLdcTTP19sgPtrrHhcaU9KQhaPFmCV1zuQQvHccoymXeFVzqZFlNzS2WYfL7PnV6q4334\nQXqukuPI/ahUhyN38s5tUsmXuYtLdbzNH5buG5o7sVSHS4K0lYfabZSOw5uWWyrD5PN9rr3PK11M\n+LCUHd7HpTrcbnB4MY2s08SWCaXs8Er5D9auk57yULsTv6V1HCeGsnPlVlzc4TmOE0PmbUtJb5TW\ncZzdRvJR2tIbr7RIm8NzHMf5peH18BzHKTe4w3Mcp9yQdocnIm1EZG7UN9JOF0jteIuivpbTRSTl\nGypE5BkRWSEiswpp1UTkbRGZLyLjo8qt6RwvR0SWRHOcXqhlXbJj1RWR90Vktoh8LiLXRnpa5lfE\neOma3z4iMinqh5onIgMiPV3zixsvLfNzdk5a1/BEJAuYB5xOovb8FKCTqtptq1Iz5ldA06Lq2id5\n/JOB9cBzqtoo0gYBq1R1UOTUf62qt6RxvDuBdapqF63b9bEOAA5Q1Rkish8wjUQ7vEtIw/yKGO9C\n0jC/aMx9VfWHqP/BROAGEq3+0vX5WeOdRprm5xRNuq/wjgO+UNVFqroJ+BeJGvbpJm3RIlWdAHy7\ng1y4N+ZQCnpmpms8SMMcVXW5qs6Ifl4PzCHRBi8t8ytiPEjTZ6iqP0Q/7g1kkXhv0/n5WeNBGr+j\nTjzpdnjbekZGLKHgC50uFHhHRKaKSExZgpRTU1XzM+JXADVLYcxrRGSmiDydylvofESkHtAEmEQp\nzK/QeJ9EUlrmJyJ7iMgMEvN4X1Vnk8b5xYwHaf78HJt0O7zdseflJFVtApwFXBXdEpYamlgjSPe8\nBwMHA42BrwG7Rv4uEt1evgz03qFLfFrmF433UjTeetI4P1XdqqqNSbT6O0VETt3h+ZTOzxgvmzR/\nfk486XZ4S2G7Vul12b6DeMqJWrqhqt8AY0jcVqebFdF6FCJyILAynYOp6kqNAJ4ihXMUkb1IOLth\nqppfWytt8ys03vD88dI5v3xUdS2JloBNKYXPr9B4zUpjfo5Nuh3eVKCBiNQTkb2Bi0j0kUwLIrKv\niFSOfq4EtAZmFf2qlFC4N2Y3MIrwpZDoP2U+HUnRHKOen08DeapauBBdWuYXN14a57d//u2jiFQk\n0dl+OumbnzlevnONSNn8nGKgqml9kLi1nAd8Adya5rEOBmZEj8/TMR4wElgG/ExiffISoBrwDjAf\nGA9UTeN4lwLPAZ8BM0n856yZorFaAFuj92969GiTrvnFjHdWGufXCPg0Gu8z4MZIT9f84sZLy/z8\nsfOHp5Y5jlNu8EwLx3HKDe7wHMcpN7jDcxyn3OAOz3GccoM7PMdxyg3u8BzHKTe4w3Mcp9zgDs9x\nnHLD/wPbeH+WuWtK7AAAAABJRU5ErkJggg==\n",
"text/plain": [
"<matplotlib.figure.Figure at 0x11362ff10>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"imshow(fit, interpolation='nearest')\n",
"title(\"Model\")\n",
"colorbar()\n",
"figure()\n",
"imshow(fit-stamp, interpolation='nearest')\n",
"title(\"Residuals\")\n",
"colorbar()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The results object also returned by py3shape.analyze contains the parameters"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"e1 = -0.0991924450855 e2 = -0.108925378414 r = 2.44430739985\n"
]
}
],
"source": [
"p=results.get_params()\n",
"print \"e1 = {} e2 = {} r = {}\".format(p.e1, p.e2, p.radius)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"You can view the full set of results like this:"
]
},
{
"cell_type": "code",
"execution_count": 12,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"x0 = 20.7746418381\n",
"y0 = 21.8365569328\n",
"e1 = -0.0991924450855\n",
"e2 = -0.108925378414\n",
"radius = 2.44430739985\n",
"radius_ratio = 1.0\n",
"bulge_A = 11.9579909796\n",
"disc_A = 0.0184697118497\n",
"bulge_index = 4.0\n",
"disc_index = 1.0\n",
"delta_e_bulge = 0.0\n",
"delta_theta_bulge = 0.0\n",
"identifier = 2016712161950\n",
"time = 0.886273\n",
"bulge_flux = 0.309212772774\n",
"disc_flux = 0.101454500524\n",
"flux_ratio = 0.752952068205\n",
"snr = 39.5264701124\n",
"old_snr = 39.5264703261\n",
"min_residuals = -0.00747191417352\n",
"max_residuals = 0.00883914731643\n",
"model_min = -0.000561836218109\n",
"model_max = 0.0273510455046\n",
"likelihood = -884.699522726\n",
"levmar_start_error = 95562.3962571\n",
"levmar_end_error = 1769.39904545\n",
"levmar_resid_grad = 0.000871920692817\n",
"levmar_vector_diff = 2.53512757471e-10\n",
"levmar_error_diff = 2.50534210886e-08\n",
"levmar_comp_grad = 9.56402914534e-10\n",
"levmar_iterations = 9\n",
"levmar_reason = 8\n",
"levmar_like_evals = 88\n",
"levmar_grad_evals = 9\n",
"levmar_sys_evals = 9\n",
"mean_flux = 0.0\n",
"number_varied_params = 7\n",
"covmat_0_0 = 0.0060332961366\n",
"covmat_0_1 = -0.00043210424842\n",
"covmat_0_2 = -1.24030315327e-05\n",
"covmat_0_3 = 1.25152574387e-07\n",
"covmat_0_4 = 0.000458099518285\n",
"covmat_0_5 = -0.0109263163141\n",
"covmat_0_6 = 1.09982900587e-05\n",
"covmat_1_0 = -0.00043210424842\n",
"covmat_1_1 = 0.00688973918684\n",
"covmat_1_2 = -6.76870715e-06\n",
"covmat_1_3 = 9.59872469977e-07\n",
"covmat_1_4 = 0.000127054462571\n",
"covmat_1_5 = -0.00467795980824\n",
"covmat_1_6 = 8.82862240801e-06\n",
"covmat_2_0 = -1.24030315327e-05\n",
"covmat_2_1 = -6.76870715e-06\n",
"covmat_2_2 = 0.00439693514216\n",
"covmat_2_3 = 7.76674301329e-05\n",
"covmat_2_4 = 0.000612640788844\n",
"covmat_2_5 = -0.0228869149194\n",
"covmat_2_6 = 6.02343539811e-05\n",
"covmat_3_0 = 1.25152574382e-07\n",
"covmat_3_1 = 9.59872469974e-07\n",
"covmat_3_2 = 7.76674301329e-05\n",
"covmat_3_3 = 0.00429349112512\n",
"covmat_3_4 = 0.00109308384598\n",
"covmat_3_5 = -0.0277226530214\n",
"covmat_3_6 = 5.98431913822e-05\n",
"covmat_4_0 = 0.000458099518285\n",
"covmat_4_1 = 0.000127054462571\n",
"covmat_4_2 = 0.000612640788844\n",
"covmat_4_3 = 0.00109308384598\n",
"covmat_4_4 = 0.119376424324\n",
"covmat_4_5 = -0.0962248542357\n",
"covmat_4_6 = -0.00477513131199\n",
"covmat_5_0 = -0.0109263163141\n",
"covmat_5_1 = -0.00467795980824\n",
"covmat_5_2 = -0.0228869149194\n",
"covmat_5_3 = -0.0277226530214\n",
"covmat_5_4 = -0.0962248542356\n",
"covmat_5_5 = 5.84685284827\n",
"covmat_5_6 = -0.020686496286\n",
"covmat_6_0 = 1.09982900587e-05\n",
"covmat_6_1 = 8.828622408e-06\n",
"covmat_6_2 = 6.02343539811e-05\n",
"covmat_6_3 = 5.98431913822e-05\n",
"covmat_6_4 = -0.00477513131199\n",
"covmat_6_5 = -0.020686496286\n",
"covmat_6_6 = 0.000298968753642\n"
]
}
],
"source": [
"for key, value in results.as_dict(0,results.number_varied_params)[0].items():\n",
" print \"{} = {}\".format(key, value)"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "IPython (Python 2)",
"language": "python",
"name": "python2"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 2
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython2",
"version": "2.7.10"
}
},
"nbformat": 4,
"nbformat_minor": 0
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment