{ "cells": [ { "cell_type": "markdown", "id": "dd69f0b0", "metadata": {}, "source": [ "# Models with Multiple Source Populations" ] }, { "cell_type": "markdown", "id": "cc4eb8de", "metadata": {}, "source": [ "*ARES* can handle an arbitrary number of source populations. To\n", "access this functionality, create a dictionary representing each source\n", "population of interest. Below, we'll create a population representative of PopII stars and another representative of PopIII stars.\n", "\n", "Before we start, it is important to note that in *ARES*, source populations are identified by their spectra over some contiguous interval in photon energy. This can be somewhat counterintuitive. For example, though UV emission from stars and X-ray emission from their compact remnants, e.g., X-ray binary systems, are both natural byproducts of star formation, we treat them as separate source populations in *ARES* even though the emission from each type of source is related to the same rate of star formation. However, because stars and XRBs have very different spectra, whose normalizations are parameterized differently, it is more convenient in the code to keep them separate. Because of this, what you might think of as a single source population (stars and their remnants) actually constitutes *two* source populations in *ARES*. \n", "\n", "Let's start with a PopII source population, and a few standard imports:" ] }, { "cell_type": "code", "execution_count": 1, "id": "b2675b7f", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Populating the interactive namespace from numpy and matplotlib\n" ] } ], "source": [ "%pylab inline\n", "import ares\n", "import numpy as np\n", "import matplotlib.pyplot as pl" ] }, { "cell_type": "code", "execution_count": 2, "id": "1a1c05d5", "metadata": {}, "outputs": [], "source": [ "pars = \\\n", "{\n", " 'problem_type': 100, # Blank slate global 21-cm signal\n", "\n", "\n", " # Setup star formation\n", " 'pop_Tmin{0}': 1e4, # atomic cooling halos\n", " 'pop_fstar{0}': 1e-1, # 10% star formation efficiency\n", " \n", " # Setup UV emission\n", " 'pop_sed_model{0}': True,\n", " 'pop_sed{0}': 'bb', # PopII stars -> 10^4 K blackbodies\n", " 'pop_temperature{0}': 1e4,\n", " 'pop_rad_yield{0}': 1e42,\n", " 'pop_fesc{0}': 0.2,\n", " 'pop_Emin{0}': 10.19, \n", " 'pop_Emax{0}': 24.6,\n", " 'pop_EminNorm{0}': 13.6,\n", " 'pop_EmaxNorm{0}': 24.6,\n", " 'pop_lya_src{0}': True,\n", " 'pop_ion_src_cgm{0}': True,\n", " 'pop_heat_src_igm{0}': False,\n", " \n", " # Setup X-ray emission\n", " 'pop_sed{1}': 'pl',\n", " 'pop_alpha{1}': -1.5, \n", " 'pop_rad_yield{1}': 2.6e38,\n", " 'pop_Emin{1}': 2e2, \n", " 'pop_Emax{1}': 3e4,\n", " 'pop_EminNorm{1}': 5e2,\n", " 'pop_EmaxNorm{1}': 8e3,\n", " \n", " 'pop_lya_src{1}': False,\n", " 'pop_ion_src_cgm{1}': False,\n", " 'pop_heat_src_igm{1}': True,\n", " \n", " 'pop_sfr_model{1}': 'link:sfrd:0',\n", "}" ] }, { "cell_type": "markdown", "id": "02ed2268", "metadata": {}, "source": [ "**NOTE:** See [problem_types](../problem_types.html) for more information about why we chose ``problem_type=100`` here. " ] }, { "cell_type": "markdown", "id": "4c17ae82", "metadata": {}, "source": [ "We might as well go ahead and run this to establish a baseline:" ] }, { "cell_type": "code", "execution_count": 3, "id": "440106f7", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "# Loaded $ARES/input/inits/inits_planck_TTTEEE_lowl_lowE_best.txt.\n", "\n", "############################################################################\n", "## ARES Simulation: Overview ##\n", "############################################################################\n", "## ---------------------------------------------------------------------- ##\n", "## Source Populations ##\n", "## ---------------------------------------------------------------------- ##\n", "## sfrd sed radio O/IR Lya LW LyC Xray RTE ##\n", "## pop #0 : fcoll yes - - x x x - ##\n", "## pop #1 : sfrd->0 yes - - - - - x ##\n", "## ---------------------------------------------------------------------- ##\n", "## Physics ##\n", "## ---------------------------------------------------------------------- ##\n", "## cgm_initial_temperature : [10000.0] ##\n", "## clumping_factor : 1 ##\n", "## secondary_ionization : 1 ##\n", "## approx_Salpha : 1 ##\n", "## include_He : False ##\n", "## feedback_LW : False ##\n", "############################################################################\n", "# Loaded $ARES/input/hmf/hmf_ST_planck_TTTEEE_lowl_lowE_best_logM_1400_4-18_z_1201_0-60.hdf5.\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ "gs-21cm: 100% |#############################################| Time: 0:00:04 \n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAEsCAYAAADTvkjJAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAAAwyElEQVR4nO3deXwV9fX/8ddJAgECyC5rjQqYsEZEcKkComUVREVxpbVI4Se2FpWKSmut2NZufsWtKMpWxKggyhIFEVCRTSAssm8l7AoCISxZzu+PO6FXyJ5779w79zwfj3nk3pm5M2/Hy8knn5n5jKgqxhhjIl+M2wGMMcYEhhV0Y4zxCCvoxhjjEVbQjTHGI6ygG2OMR1hBN8YYj7CCbowxHmEF3RhjPMIKujHGeIQVdGMCRESGiYgWMM12O5uJDmK3/hsTGCJSFajqN+t24K/ALao6151UJppYQTcmCETkXuAV4FZV/cztPCY6WJeLMQEmIoOAMUAvK+YmlKygGxNAIjIMXzdLN1X90u08JrrEuR3AGK8QkceA3wE3quoqt/OY6GN96MYEgIiMBJ4AbgXW+y3KVNVMd1KZaGMF3ZhyEhEBfgCqF7D4EVX9v9AmMtHKCroxxniEnRQ1xhiPsIJujDEeYQXdGGM8wgp6gInIb0VkvYisE5F3RKSSiNQSkbkissX5WdOlbJVEZJmIpDsZ/+jMD0k+EWkiIp+LyAZn/78J1f5F5C0ROSgi6/zmPSMie0RktTP1DME+/yYiG0VkjYhMF5EaIdhnWHz/ClPQvxm3M4WaiNQQkfed78YGEbm6LNuxgh5AItII+DXQXlVbAbHAAHyXs32mqs2Az5z3bjgN3KCqbYEUoLuIXBXCfDnAo6qaDFwFPCQiLUK0//FA9wLm/0tVU5wp0INoFbTPuUArVW0DbAZGhmCf4fL9O08R/2aizf8BaaqaBLQFNpRlI1bQAy8OqCwicUAVYC/QF5jgLJ8A3OJGMPXJvya6gjMpIcqnqvtUdaXz+ji+L22jUOxfVRcBhwO93dLuU1U/VdUc5+0SoHGw90mYfP+KUNC/maghItWB64FxAKp6RlV/KMu2rKAHkKruAf4O/BfYBxxV1U+BC1V1n7POPqCeWxlFJFZEVgMHgbmqutSNfCKSCFwOuLJ/P8Oc7o+3XOiKeACYE4L9hM3371xF/JuJJpcAh4C3RWSViLwpIgll2ZAV9AByCkJf4GKgIZDgjLoXNlQ1V1VT8LUMO4hIq1BncIaZ/QDfTTfHQr1/P68Bl+LrftoH/CNUOxaRp/B1Qf0nVPsMR5HwbyYE4oB2wGuqejlwgjJ2i1lBD6wbgR2qekhVs4FpwDXAARFpAOD8POhiRgCcP+kW4OtvDVk+EamAr5j/R1WnObNdOT6qesD5BZcHvAF0CMV+RWQg0Bu4R0NzZ1/Yff/8FPZvJppkABnOX8sA7+Mr8KVmBT2w/gtcJSJVnNvBu+LrJ/4IGOisMxCY4UY4Eambf1WFiFTG949pY6jyOcdkHLBBVf/pt8iV45Nf5Bz9gHWFrRvAfXbHN4BXH1XNCvb+HGHx/StEYf9mooaq7gd2i8hlzqyuwLdl2Zbd+h9gzqWAd+L7c3oVMAjfU2xSgZ/g+wL3V9WQnqBzsrXBd1IsFt8v81RVfVZEaocin4j8FPgCWAvkObOfxNePHtT9i8g7QGegDnAA+IPzPgXfieGdwK/y+5qDuM+RQDzwvbPaElUdEuR9fkgYfP8KU9C/GVU97W6q0BKRFOBNoCKwHfiFqh4p9XasoBtjjDdYl4sxxniEFXRjjPEIK+jGGOMRVtCNMcYjIrKgi8hOEVnrDKi0wpkXlgMQichgtzMUx82M0bbvaNlnaUVCxmALxDGIyILu6OIMqNTeeR+uAxCV+X+SiNwcyPWK4OY/pmjbd7Tss7QiIWOwRXVBP1e4D0BUFiUt1OUt6MYYD4jI69BFZAdwBN8NIf9W1bEi8oOq1vBb54iqntft4vxZk/+b8IoqVaoENWtOTg5xcXFl+mxubi6xsbEBW68w5clYXtG272jZZ2lFQsZgK+kxyMrKUlUtsDEeqQW9oaruFZF6+MaXfhj4qCQF3V9CQoKeOHEiuGGNMSaARCRLVQscjTEiu1xUda/z8yAwHd+gSuE8AJExxgRdxBV0EUkQkWr5r4Gf4RtUKZwHIDLGmKCLxE6rC4HpvoHZiAOmqGqaiCwHUkXklzgDELmY0RhjQi4i+9ADxfrQjTGRxnN96MYYY85nBd0YYzzCCroxxniEFXRjjPEIK+jGGOMRVtCNMcYjrKAbY4xHWEE3xhiPsIJujDEeYQXdGGM8wgq6McZ4hBV0Y4zxCCvoxhjjEVbQjTHGI6ygG2OMR1hBN8YYj7CCbowxHmEF3RhjPMIKujHGeIQVdGOM8Qgr6MYY4xFW0I0xxiOsoBtjjEdYQTfGGI+wgm6MMR5hBd0YYzwizu0AxpjopaqcOHGCI0eOkJmZSVZWVoHTyZMnOX36NDk5OT+asrOzz5unqgCIyNn95L8u7GdMTAyxsbHExcUF7Gf+VKFCBSpUqHD29bk/i1rm/99QEp4q6CLSHfg/IBZ4U1X/4nIkY6JOdnY2+/btY+/evezZs4e9e/eenb777juOHDnC4cOHOXLkCEeOHCE7O7tM+8kvfP5TbGwsMTG+jof8wl7cT4Dc3NyzU05Oztmf/uu4If+Xg3+RL4q4HThQRCQW2AzcBGQAy4G7VPXbwj6TkJCgJ06cCFFCY7xDVcnIyGDNmjVs2rSJrVu3sm3bNrZu3cquXbvIzc390fpxcXE0aNCAunXrUqtWLWrWrHl2yn9frVo1qlSpcnZKSEg4+7py5crEx8f/qHCHQl5eXoGFvrifBf0V4f+zoHklXfb6669nqWpCQXm9VNCvBp5R1W7O+5EAqvrnwj5TuXJl3bp1K7Vq1aJy5cohSmpMZFFVdu3axVdffcXy5ctJT08nPT2dI0eOnF2nRo0aNG3alKZNm3LppZeSmJhIo0aNaNiwIQ0bNqR27dpnW86mfESk0ILupS6XRsBuv/cZQMeiPnDq1CkaN25MTEwMTzzxBM8991yp+6yM8aJt27aRlpbGwoUL+eqrr9i7dy8AVapUoXXr1vTv3582bdrQtm1bWrRoQa1atVxObMBbBb2gSnzenx8iMhgYDL4/A19++WXmz5/P888/T+3atRk+fHiwcxoTdnJycliwYAEzZ85kzpw5bN68GYAmTZrQqVMnrr32Wn7605/SqlWrkHV3mNKL6i6X/D70vLw8br/9dmbOnMnKlStp1apViFIb4568vDwWL17MO++8w3vvvcehQ4eoVKkSnTt3pmfPnvTo0YOmTZu6HdOco6guFy8V9Dh8J0W7AnvwnRS9W1XXF/YZ/5Oihw4domXLliQlJbFw4ULrejGe9f333/P222/z+uuvs23bNipVqsTNN9/MXXfdRbdu3ahSpYrbEU0RiironjlLoao5wDDgE2ADkFpUMT9X3bp1+eMf/8gXX3zB7NmzgxXTGNds2LCBX/ziFzRq1IjHH3+cBg0aMHHiRA4ePEhqair9+vWzYh7hPNNCL4tzL1vMzs4mOTmZhIQEVq9e7YlW+nfffcf8+fNZtWoVBw4coFq1anTs2JG+ffuSkFDgL3njMatWreL555/ngw8+oHLlygwcOJChQ4fSunVrt6OZMiiqhY6qRu1UpUoVPdf48eMV0NmzZ5+3LFLk5OTohx9+qD179lQRUUArVKigjRo10qpVqyqg9evX148//tjtqCaItm/frnfccYcCWr16dX3qqaf04MGDbscy5QSc0EJqmutF1c2poIJ++vRpbdy4sXbq1KlUBzkc5Obm6tSpUzUpKUkBbdiwoT755JP69ddfa3Z29tl1Pv/8c7388stVRHTy5MkupzaBdvToUR0xYoRWrFhRK1eurKNGjdIjR464HcsEiBX0UhR0VdV//OMfCuiSJUtKfpRdNmfOHG3durUC2rJlS01NTT1bxAuSlZWlnTt31goVKujKlStDmNQE00cffaSNGjVSEdGBAwdqRkaG25FMgFlBL2VBP3bsmNaoUUNvvfXWkh9ll+zatUv79eungDZt2lSnTJmiOTk5Jfrs999/r/Xr19eUlJQSf8aEp0OHDumAAQMU0NatW+vSpUvdjmSCxAp6KQu6quqTTz6pIqKbN28u2VF2wYQJEzQhIUErV66szz//vJ46darU25g6daoCOmXKlCAkNKHw+eefa8OGDbVixYr67LPP6unTp92OZILICnoZCvq+ffu0YsWKOmTIkJId5RA6fvy4Dhw4UAHt1KmT7ty5s8zbys3N1ZYtW2pycrLm5uYGMKUJtpycHH3mmWc0JiZGmzdvrqtWrXI7kgkBK+hlKOiqqoMGDdL4+Hg9cOBA8Uc5RNLT0zUpKUlFRP/whz8EpKtk8uTJCujcuXMDkNCEwtGjR7VHjx4K6H333afHjx93O5IJESvoZSzoGzZsUEBHjRpV/FEOgbfeeksrVaqk9evX188++yxg2z158qTWrFlTBwwYELBtmuDZvn27tmjRQuPi4vT11193O44JMSvoZSzoqqp9+vTRWrVqaWZmZrHrBktOTo7+9re/VUC7du2q+/fvD/g+hg0bpvHx8Xr06NGAb9sEzooVK7ROnTpas2bNgP5SN5GjqILumVv/g2XEiBEcPnyYt99+25X9Z2Zm0qdPH/71r3/x8MMPk5aWxoUXXhjw/dx5552cPn2atLS0gG/bBMaiRYvo0qULCQkJLFmyhBtuuMHtSCbM2K3/JXhi0TXXXMP+/fvZvHkzcXGhG3H4hx9+oGfPnixdupRXXnmFIUOGBG1fubm5NGjQgBtvvJEpU6YEbT+mbNLS0ujXrx+JiYnMnTuXxo0bux3JuCQqBucKpscff5wdO3Ywbdq0kO3z4MGDdOnShRUrVpCamhrUYg6+Zxf27t2b2bNnc+bMmaDuy5TOvHnz6Nu3L8nJySxatMiKuSmUFfQS6NOnD82aNeOFF14gFH/RZGRkcP3117Np0yY+/vhjbrvttqDvE6B3794cPXqU5cuXh2R/pniLFy+mb9++NG/enHnz5lG3bl23I5kwZgW9BGJjY3niiSf45ptv+PDDD4O6r23btnHdddexd+9ePvnkE7p16xbU/fnr1KkTAAsWLAjZPk3hVq9eTc+ePWnYsCFz5861x7yZYlkfegn60MH3iK784UbXrl0blL709evXc9NNN3HmzBnS0tJo3759wPdRnLZt21KvXj3mzp0b8n2b/9mzZw8dOnQgJiaGL7/8kosuusjtSCZMWB96AMTFxfH888+zceNGxo8fH/Dtf/PNN2dbyAsXLnSlmAN07tyZr776itOnT7uyfwNZWVn07duXY8eOMWvWLCvmpsSsoJfCLbfcwlVXXcXvf/97jh49GrDtfvHFF9xwww1UrVqVL774gpYtWwZs26XVpUsXTp48af3oLsnLy+P+++9n5cqVvPPOO7Rp08btSCaCWEEvBRFhzJgxHDhwgBEjRgRkmzNnzqRbt240aNCAL7/8kksvvTQg2y2ra665BoAlS5a4miNa5T9Z6O9//zu9e/d2O46JNIXdcRQNU0nuFC3IY489pkC579R7/fXXNSYmRtu3bx9W48UkJiZq//793Y4RdebPn68xMTF69913a15enttxTJiiiDtF7aRoCU+K+svKyqJdu3b88MMPrFy5koYNG5bq86rKqFGjGD16ND179iQ1NTWsnu85YMAAlixZws6dO92OEjUOHDhASkoKF1xwAStWrKBq1apuRzJhyk6KBliVKlX44IMPyMzM5PbbbycrK6vEnz169Ch33nkno0ePZtCgQcyYMSOsijlAhw4d2LVrFwcOHHA7SlTIzc3lnnvu4ejRo7z33ntWzE2ZWUEvo5YtWzJ+/HiWLFlC7969yczMLPYzy5Yt4/LLL2fatGn89a9/ZezYsSEdSqCkOnToAMDSpUtdThIdXnzxRT777DPGjBlz9tJYY8qksL6YaJjK2ofub/LkyRoTE6OXXHKJTp8+/byHROTl5enq1at12LBhGhsbqxdddJEuXry43PsNpszMTBURfeaZZ9yO4nnffvutxsfHa9++fa3f3JQI1odesLL2oZ9r0aJFDBo0iC1btlCvXj2uuuoqatSowcGDB1m3bh0ZGRnExMTwq1/9itGjR1OzZs0ApA+uyy67jJYtW4Z0/Jpok5OTwzXXXMP27dtZv359UEbRNN5TVB96+P29H4Guv/56vv32W6ZNm8ZHH31Eeno6x44do3bt2lx77bXcdNNN9OrVi/r167sdtcTatGnDqlWr3I7haS+88ALLly/n3XfftWJuAsJa6AFooXvRc889x6hRozh27BjVqlVzO47nbNq0iTZt2tC3b19SU1PdjmMiiF3lYkot/w7FtWvXupzEe1SVYcOGUblyZV566SW34xgPsYJuCtS2bVsA1qxZ43IS73nvvfeYN28eo0ePjqhuOBP+Iqqgi8gzIrJHRFY7U0+/ZSNFZKuIbBKR0I0561E/+clPuOCCC0hPT3c7iqccO3aMRx55hHbt2gX9oSUm+kTiSdF/qerf/WeISAtgANASaAjME5HmqprrRkAvEBHatGljLfQAe+aZZ9i/fz8ffvghsbGxbscxHhNRLfQi9AWmquppVd0BbAU6uJwp4rVp04a1a9cSzSfOA2nLli2MGTOGQYMGnb15y5hAisSCPkxE1ojIWyKSf0F3I2C33zoZzrzziMhgEVkhIitycnKCnTWitWjRguPHj7Nnzx63o3jCyJEjiY+P509/+pPbUYxHhV1BF5F5IrKugKkv8BpwKZAC7AP+kf+xAjZVYLNSVceqantVbR+Ot92Hk+TkZAA2bNjgcpLIt3jxYj744ANGjBhh15yboAm7iqaqN5ZkPRF5A5jpvM0AmvgtbgzsDXC0qONf0G+66SaX00QuVeWxxx6jQYMGPProo27HMR4Wdi30oohIA7+3/YB1zuuPgAEiEi8iFwPNgGWhzuc1F154ITVq1LAWejlNmzaNr7/+mmeffTbsRtY03hJRd4qKyCR83S0K7AR+par7nGVPAQ8AOcAjqjqnuO3ZnaLFu/rqq4mPj2fBggVuR4lIubm5tGrVitjYWNLT0+3KFlNunhnLRVXvK2LZaGB0CONEheTkZGbNmuV2jIj17rvvsnHjRt5//30r5iboSt3lIiKXiUgPEblVRK4TERuN38OSk5M5ePAghw8fdjtKxMnNzeXZZ5+ldevW9OvXz+04JgqUqIUuIonA/wPuBS7kx1eV5IjIl8DrwHsaSX04plj+J0avvfZal9NElnfffZdNmzbx/vvvExMTUaerTIQq9lsmIn/Dd/LxMuBJoBVwARAPNAB6AouBvwCrRaRd0NKakLNLF8vGWufGDSVpoVcDmqtqQZcBHnCmecDTItIfSAZWBi6icVNiYiLx8fFs3LjR7SgRxVrnxg0RdZVLoNlVLiXTtm1bGjdubCdHS0hVzw4/nJ6ebgXdBFS5x0MXkfgSrJNSylwmQiQnJ1uXSymkpaWxbt06RowYYcXchFRJv21TilrojHb4SfnjmHCUnJzMzp07ycrKcjtKRHjhhRdo3LgxAwYMcDuKiTIlLejXiMiYghaISFN8feh2Z6ZHJSUloaps2bLF7Shhb9myZSxYsIDf/va3VKhQwe04JsqUtKD3AgaKyEj/mSLyE+AzYANwW4CzmTCRf6WLnRgt3t/+9jcuuOACHnzwQbejmChUooKuqiuB24E/iMh9cHZclc+APcDNqnomaCmNq5o1a4aIWD96MbZu3cq0adMYOnSoPVjbuKLEt/6r6qci8iDwhojk4bsm/RjQQ1Wtc9XDKleuTGJiorXQi/Hiiy8SFxfHr3/9a7ejmChVqrFcVHWSiNQHJuLrZrlJVY8GJZkJK8nJyVbQi3Ds2DEmTJjAnXfeSYMGDYr/gDFBUNJb/z89Z1a2M00V+d8oAKr6s8BFM+EkKSmJ+fPnk5uba4NMFWDChAlkZmby8MMPux3FRLGSttDPfQbZO4EOYsJbUlISp06d4r///S8XX3yx23HCSl5eHi+//DIdO3bkyiuvdDuOiWIlKuiq+otgBzHhzf9KFyvoPzZv3jw2b97MpEmT3I5iopzdxmZKJCkpCbBBugoyZswY6tWrR//+/d2OYqJcqR9wISLdgK5APc75haCq9wcolwkzderUoXbt2nZi9Bzbt29n1qxZPP3008THFztChjFBVaqCLiLP4btccQ2wH9+j4EyUsCtdzvfqq68SGxvLkCFD3I5iTKlb6IOBn6vqxGCEMeEtKSmJDz/80O0YYePUqVO8/fbb9OvXj4YNG7odx5hS96Hn4XuYhYlCSUlJfPfdd3z33XduRwkL06dP5/DhwwwePNjtKMYApS/orwKDghHEhL/8K102bdrkcpLw8Oabb3LxxRdzww03uB3FGKD0XS5/AmaKSDq+fvRs/4Wq+kCggpnwk3+ly8aNG6P++aLbtm1j/vz5PPfcczbmuQkbpS3ozwI9gG/xPU/UTopGkYsuuoj4+Hi7dBFf6zw2NpZf/MJu0TDho7QFfRjwgKqOD0IWE+ZiY2Np3rx51F/pkp2dzfjx4+nVq5edDDVhpbR/K54BvgxGEBMZ7NJFmDVrFvv372fQIDudZMJLaQv6WOCXwQhiIkNSUhI7duzg1KlTbkdxzRtvvEHDhg3p0aOH21GM+ZHSdrk0AG5z7hZN5/yTonb9lsclJSWRl5fHli1baN26tdtxQm737t2kpaXx5JNPEhdX6hutjQmq0rbQLwVWA0eBRKCZ39Q0UKFEpL+IrBeRPBFpf86ykSKyVUQ2Ob9Y8udfISJrnWUvif+4viZgov1xdG+//Taqyi9/aX+omvBT2gdcdAlWkHOsA24F/u0/U0RaAAOAlkBDYJ6INFfVXOA1fHeyLgFmA92BOSHKGzWaN28OROcgXarK+PHj6dq1K4mJiW7HMeY8xbbQReSKkm5MRCqJSHL5IoGqblDVgu5e6QtMVdXTqroD2Ap0cJ5vWl1Vv1ZVxfdEpVvKm8Ocr0qVKlx00UVR2UL/8ssv2bFjBwMHDnQ7ijEFKkmXywwRmS4i3USkwPVFpJGIjAS2AMG846QRsNvvfYYzr5Hz+tz55xGRwSKyQkRW5OTkBC2ol0XrlS4TJ06katWq9OvXz+0oxhSoJF0ulwFPAJOBSiKyCt8TjE4BtfB1f1wMLADuUtUSXdYoIvOA+gUsekpVZxT2sQLmaRHzz5+pOhbf1TokJCTYjVFlkJSUxMKFC8nLy4uauyRPnjxJamoqt99+OwkJCW7HMaZAxRZ0VT0BjHKGzu0BXI+vgNcEDgGvAHNUtVRNNlW9sfRxyQCa+L1vDOx15jcuYL4JgqSkJE6ePMnu3bu56KKL3I4TEjNmzODYsWPcf78N+W/CV4lPiqrqaeBDZ3LLR8AUEfknvpOizYBlqporIsdF5CpgKXA/MMbFnJ7mf6VLtBT0iRMn0qRJEzp16uR2FGMKFZZ/L4tIPxHJAK4GZonIJwCquh5IxTeWTBrwkHOFC8BQ4E18J0q3YVe4BE20PY5u3759fPLJJ9x3331R08VkIlNY3hmhqtOB6YUsGw2MLmD+CqBVkKMZoG7dutSsWTNqToxOmTKFvLw8624xYc+aG6bURCRqrnRRVSZMmEDHjh257LLL3I5jTJGsoJsySUpKioqCnp6eztq1a611biKCFXRTJklJSRw4cIAjR464HSWoJk6cSIUKFRgwYIDbUYwplhV0Uyb5V7p4+cRodnY2//nPf7j55pupVauW23GMKVa5CrqI/ExEVorILhGZISLtAhXMhLcWLVoAsH79epeTBM+nn37KwYMHrbvFRIzyttBfA36D7+qSfwH/FJH7yp3KhL3ExESqV6/O6tWr3Y4SNBMnTqROnTo27rmJGOUt6AdV9QtVPa6qC4CewO/KH8uEu5iYGNq2bevZgn7kyBFmzJjBXXfdRcWKFd2OY0yJlKmgi8i7IjIC+EpE/iQiFZxFucDpgKUzYS0lJYX09HTy8vLcjhJw7733HqdPn7buFhNRytpCHwNkATXwje+yTUQWAJvxjUVuokBKSgonTpxg+/btbkcJuAkTJtCiRQuuuKLEo0cb47oS3SkqIm8Av1HVLABnRMUv/ZbHAsnA5UCbIOQ0Yaht27YArF69mqZNA/bAKtdt3bqVxYsX8+c//xl78JWJJCVtoT8AVC1soarmquo6VZ2kqo8HJpoJdy1btiQ2NtZz/eiTJk1CRLj33nvdjmJMqZS0oFszxZynUqVKJCcne6qgqyqTJk2ia9euNG7cuPgPGBNGStOHbg+DMOfJPzHqFV999RU7duywk6EmIpWmoP9bRB4XkS4iUj1oiUxESUlJISMjg0OHDrkdJSAmTpxIlSpV7DFzJiKVpqA3BJ4GPgOOiMhmEZkiIsNFpJOIVAtORBPO2rdvD8CyZctcTlJ++Y+Zu+2226hatdBTRsaErdIU9D74LlNsge+JQLPwPQ7uj8DngLdHaTIFuuKKK4iJiWHp0qVuRym3jz/+mKNHj1p3i4lYJX3AhQKoqgIbnek/AOK7rqsFYOO4RKGqVavSqlUrT7TQJ02aRKNGjejSpYvbUYwpk3Jf5aI+61V1UoAymQjTsWNHli1bhu/3fWQ6ePAgc+bM4d577yU2NtbtOMaUSUkLelfghyDmMBGsY8eOHDlyhC1btrgdpczeeecdcnNzue8+G1vORK6SFvRhQF8AEUkSkU9EZJ2ITBWRnsGLZyJBx44dASK6H33ixIm0a9eOli1buh3FmDIraUHvBKxxXk8BDgNfAJcAM0XkLbF7pKNWcnIyVatWjdiCvn79elauXGknQ03EK+lJ0arAMRFpA7ypqq/mLxCR64BUYDDw78BHNOEuNjaWK6+8kiVLlrgdpUwmTZpEbGwsd911l9tRjCmXkrbQDwANgC74ivdZqvoF8Ajwq4AmMxHlmmuuYfXq1Rw/ftztKKWSm5vL5MmT6d69O/Xq1XM7jjHlUtKCPht4A3gcKKiT8RvAO8PtmVLr3Lkzubm5fPXVV25HKZW5c+eyZ88eBg4c6HYUY8qtpAX9d8BKfDcQJYnIEBHxv5XuHmBfoMOZyHH11VdToUIFFixY4HaUUhk3bhy1a9emT58+bkcxptxK1IeuqseAB/Pfi8ifge9EZBeQgK875rdBSWgiQkJCAh06dIiogn7o0CFmzJjBQw89RHx8vNtxjCm3Mj2xSFVHAinABHx96neo6ksBzGUiUOfOnVmxYkXE9KNPnjyZ7OxsfvnLX7odxZiAKPNDolV1o6o+r6rDVfWDQIYSkf4isl5E8kSkvd/8RBE5KSKrnel1v2VXiMhaEdkqIi/ZZZShF0n96KrKuHHj6NChA61atXI7jjEBUeaCHmTrgFuBRQUs26aqKc40xG/+a/gunWzmTN2DH9P4y+9Hnz9/vttRirVs2TLWr19vrXPjKWFZ0FV1g6puKun6ItIAqK6qXzsDiE0EbglWPlOwhIQErrvuOmbPDv/nhI8bN47KlSszYMAAt6MYEzBhWdCLcbGIrBKRhc5NTQCNgAy/dTKceecRkcEiskJEVuTk5AQ7a9Tp1asX69evZ9euXW5HKdSJEyeYOnUq/fv3p3p1e1aL8Q7XCrqIzHPGgzl36lvEx/YBP1HVy4HhwBTn6UkF9ZcXOPSfqo5V1faq2j4urqQ3ypqS6tWrF0BYt9Lfeecdjh8/zqBBg9yOYkxASTgPeSoiC4DHVHVFUcuBPcDnqprkzL8L6KyqRd69mpCQoCdOnAho5minqjRt2pTk5GRmzpzpdpzzqCrt2rUjJyeHNWvWYOfOTaQRkSxVTShoWUR1uYhIXRGJdV5fgu/k53ZV3QccF5GrnKtb7gdmuBg1aokIvXr1Yv78+Zw8edLtOOdZunQpq1ev5qGHHrJibjwnLAu6iPQTkQzgamCWiHziLLoeWCMi6cD7wBBVPewsGwq8CWwFtgFzQhzbOHr37s3JkyeZO3eu21HO88orr1CtWjXuuecet6MYE3Bh3eUSbNblEhxnzpyhfv369OzZk8mTJ7sd56xDhw7RuHFjBg8ezJgxY9yOY0yZeKbLxUSGihUrcuuttzJjxoyw6nYZN24cZ86cYejQoW5HMSYorKCboLjzzjvJzMwkLS3N7SgAZGdn8+qrr9KlSxdatGjhdhxjgsIKugmKLl26UKdOHaZOnep2FABSU1PZvXs3jz76qNtRjAka60O3PvSgefjhh3njjTfYu3cvtWrVci2HqnL55ZeTnZ3N2rVriYmxdoyJXNaHblwxaNAgTp8+zaRJk1zN8dlnn5Gens6jjz5qxdx4mrXQrYUeVB06dODkyZOu3sTTrVs31qxZw86dO23ccxPxrIVuXPPggw+ybt061x4gvXz5cj799FN+/etfWzE3nmctdGuhB9Xx48dp0qQJP/vZz0hNTS3+AwHWq1cvlixZwo4dO2wgLuMJ1kI3rqlWrRpDhw7lgw8+YNu2bSHd99KlS5k9ezaPPfaYFXMTFayFbi30oNu3bx+JiYkMGjSIV155JWT77dGjB8uXL2fHjh1Uq1YtZPs1JpishW5c1aBBA+677z7eeust9u7dG5J9zp8/n7S0NEaMGGHF3EQNa6FbCz0ktm/fTlJSEj//+c8ZO3ZsUPeVm5tLu3btOHbsGBs2bKBSpUpB3Z8xoWQtdOO6Sy65hKFDhzJu3Dg2bNgQ1H299dZbrFmzhhdeeMGKuYkq1kK3FnrIHDp0iEsvvZROnTrx8ccfB2Uf33//PS1atKB58+YsWrTIxjw3nmMtdBMW6taty6hRo5g5cybTpk0Lyj6GDx/O4cOHefnll62Ym6hjLXRroYdUTk4OV155Jfv37+fbb7+lZs2aAdv2nDlz6NmzJ08//TR/+tOfArZdY8JJUS10K+hW0EPum2++4aqrrqJPnz68//77AWlJ79u3j5SUFOrWrcs333xjd4Uaz7IuFxNWrrjiCp5//nmmTZsWkCcH5eTkcNddd5GZmUlqaqoVcxO1rKAbVzz66KPcfPPNDB8+nFmzZpV5O6rKQw89xMKFC3nttdfs4RUmqllBN66IiYlhypQptG3bljvuuIMFCxaUehuqyu9//3vGjh3LyJEjuf/++wMf1JgIYgXduKZq1arMmjWLxMREunfvXqorX3Jzcxk+fDjPPfccDzzwAKNHjw5iUmMigxV046r69euzaNEiUlJSuO222/jNb37D8ePHi/zM7t27uemmm3jxxRd55JFHeOONN+wSRWOwgm7CQO3atVm4cCEPP/wwL730Es2aNeMvf/kLu3btIv8qLFVlzZo1DB8+nObNm7Ns2TLGjRvHP//5T3sKkTEOu2zRLlsMK0uXLmXkyJF8/vnnANSpU4fq1atz8OBBMjMziY2N5e677+bZZ58lMTHR3bDGuMCuQy+EFfTwtWXLFtLS0li7di1ZWVnUqlWLtm3bcvPNN1OvXj234xnjGivohbCCboyJNHZjkTHGRIGwLOgi8jcR2Sgia0RkuojU8Fs2UkS2isgmEenmN/8KEVnrLHtJ7LIHY0yUCcuCDswFWqlqG2AzMBJARFoAA4CWQHfgVRGJdT7zGjAYaOZM3UMd2hhj3BSWBV1VP1XVHOftEqCx87ovMFVVT6vqDmAr0EFEGgDVVfVr9Z0UmAjcEurcxhjjprAs6Od4AJjjvG4E7PZbluHMa+S8Pne+McZEjTi3diwi84D6BSx6SlVnOOs8BeQA/8n/WAHraxHzC9rvYHxdM1SsWLGUqY0xJny5VtBV9cailovIQKA30FX/d21lBtDEb7XGwF5nfuMC5he037HAWPBdtlim8MYYE4bCsstFRLoDvwP6qGqW36KPgAEiEi8iF+M7+blMVfcBx0XkKufqlvuBGSEPbowxLnKthV6Ml4F4YK5z9eESVR2iqutFJBX4Fl9XzEOqmut8ZigwHqiMr899znlbNcYYD7M7Re1OUWNMBLE7RY0xJgpYQTfGGI+wgm6MMR5hBd0YYzzCCroxxniEFXRjjPEIK+jGGOMRVtCNMcYjrKAbY4xHWEE3xhiPsIJujDEeYQXdGGM8wgq6McZ4hBV0Y4zxCCvoxhjjEVbQjTHGI6ygG2OMR1hBN8YYj7CCbowxHmEF3RhjPMIKujHGeIQVdGOM8Qgr6MYY4xFW0I0xxiOsoBtjjEdYQTfGGI+wgm6MMR5hBd0YYzwiLAu6iPxNRDaKyBoRmS4iNZz5iSJyUkRWO9Prfp+5QkTWishWEXlJRMS1/wBjjHFBWBZ0YC7QSlXbAJuBkX7LtqlqijMN8Zv/GjAYaOZM3UOW1hhjwkBYFnRV/VRVc5y3S4DGRa0vIg2A6qr6taoqMBG4JbgpjTEmvMS5HaAEHgDe9Xt/sYisAo4BT6vqF0AjIMNvnQxn3nlEZDC+lnz++6wy5ooFciPgc+X5bByQU+xagdtfeT5rxyf89mnHp2hlPT6VC12iqq5MwDxgXQFTX791ngKmA+K8jwdqO6+vAHYD1YErgXl+n7sO+LgEGVaUI//YSPhcOfdZpuPjUlY7PuG3Tzs+QTg+RU2utdBV9cailovIQKA30FWd/3pVPQ2cdl5/IyLbgOb4WuT+3TKNgb3ByO3n4wj5XHk/G+r92fEJzmft+ITfPgMuv+UbVkSkO/BPoJOqHvKbXxc4rKq5InIJ8AXQWlUPi8hy4GFgKTAbGKOqs4vZzwpVbR+0/5AIZ8enaHZ8imbHp2jBOD7h2of+Mr7ulbnO1YdL1HdFy/XAsyKSg6/PaoiqHnY+MxQYj69/aY4zFWdsgHN7jR2fotnxKZodn6IF/PiEZQvdGGNM6YXlZYvGGGNKzwq6McZ4RNQWdBHpLiKbnKECnnA7TzgQkZ3O8AmrRWSFM6+WiMwVkS3Oz5pu5wwVEXlLRA6KyDq/eYUeDxEZ6XyfNolIN3dSh04hx+cZEdnjNzxHT79lUXN8RKSJiHwuIhtEZL2I/MaZH9zvT6Cvg4yECd+NANuAS4CKQDrQwu1cbk/ATqDOOfNeAJ5wXj8B/NXtnCE8HtcD7YB1xR0PoIXzPYoHLna+X7Fu/ze4cHyeAR4rYN2oOj5AA6Cd87oaviFMWgT7+xOtLfQOwFZV3a6qZ4CpQF+XM4WrvsAE5/UEomhIBVVdBBw+Z3Zhx6MvMFVVT6vqDmArvu+ZZxVyfAoTVcdHVfep6krn9XFgA76714P6/YnWgt4I312m+QodKiDKKPCpiHzjDJEAcKGq7gPflxSo51q68FDY8bDv1P8Mc0ZKfcuvSyFqj4+IJAKX47tHJqjfn2gt6AUNrWvXb8K1qtoO6AE8JCLXux0ogth3yuc14FIgBdgH/MOZH5XHR0SqAh8Aj6jqsaJWLWBeqY9PtBb0DKCJ3/tQDBUQ9lR1r/PzIL4xdDoAB5zRLPNHtTzoXsKwUNjxsO8UoKoHVDVXVfOAN/hft0HUHR8RqYCvmP9HVac5s4P6/YnWgr4caCYiF4tIRWAA8JHLmVwlIgkiUi3/NfAzfIOlfQQMdFYbCMxwJ2HYKOx4fAQMEJF4EbkY35j8y1zI56r8YuXoh+87BFF2fJwH7IwDNqjqP/0WBfX7E663/geVquaIyDDgE3xXvLylqutdjuW2C4HpzlALccAUVU1zxshJFZFfAv8F+ruYMaRE5B2gM1BHRDKAPwB/oYDjoarrRSQV+BbfkKgPqWpZh2ONCIUcn84ikoKvu2An8CuIyuNzLXAfsFZEVjvzniTI3x+79d8YYzwiWrtcjDHGc6ygG2OMR1hBN8YYj7CCbowxHmEF3RhjPMIKujHGeIQVdGOM8Qgr6MaUgIhUdcb5vtLtLAAi8m8R+bvbOUx4sYJuTMn8DlihqsvzZ4jIeBFREfng3JVF5BZnWc45688raOPOuveWIs+zwFARuaQUnzEeZwXdmGKISCVgKPDvAhb/F7hZRC48Z/5gYFewMqnqHuAz4P8Fax8m8lhBN1FHRN4TkU/83tcVkSwRuaKQj3QHKgOfFrBsC7AE+Lnf9n4C3AS8XcZ8nZ0W+7nTznNWnQ6UplVvPM4KuolGPxqqVFUP4Xv4wC2FrN8JWKWqOYUsHwsMckbYAxiEr/Vc1hb6YnyPMMufWuIbSvXzc9ZbClwoIsll3I/xGCvoJhqdO/Y0QCaFP43pYmBPEdt7H6iFb6TBWOABfEW+IJ1FJPPcyX8FVT2jqvtVdT/wPfAKsB0YUsB/B/iejWtMdA6fa6JeBlBVRC5Q1aPO+O8/BR4qZP3KwNHCNqaqp0RkEvAgvgcCxwEfA/cUsPpS/jcetr8thWz+NXy/fK5S1dPnLDvll88YK+gmKuW3bJvgK9R/BvYD0wpZ/xC+FnhR/g2sAn4CvK2q2f/rgfmRk6q69dyZBa0rIiOAW4GrVfW7AraVn+lQMdlMlLCCbqJRfkFv7FxXfh/wU1U9Vcj6K4FhRW1QVTc4DwO5loJb4KUiIrfguzSxu6puKmS11kAuvl8kxlhBN1FpD5AH/AbfMy9/VswTq+YA/xCRJqq6u4j1ugGVVPVwecKJSEtgMvAMsFFE6juLcp0TuPk6A18W8/BhE0XspKiJOs7VKgeAtkBX/5uFCll/A7AAX0u+qPWyylvMHVcCCfi6gvb5Tf43NQlwNwVfG2+ilD2CzpgSEJHrgKlAM1XNCoM8dwCjgBSPP5vTlIK10I0pAVX9AvgjvksYw0E88Asr5saftdCNMcYjrIVujDEeYQXdGGM8wgq6McZ4hBV0Y4zxCCvoxhjjEVbQjTHGI/4/Q4HOhf8d1VoAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "sim = ares.simulations.Global21cm(**pars)\n", "sim.run()\n", " \n", "ax, zax = sim.GlobalSignature(color='k')" ] }, { "cell_type": "markdown", "id": "f4d90c91", "metadata": {}, "source": [ "Now, let's add a PopIII-like source population. We'll assume that PopIII sources are brighter on average (in both the UV and X-ray) but live in lower mass halos. We could just copy-pase the dictionary above, change the population ID numbers and, for example, the UV and X-ray ``pop_rad_yield`` parameters. Or, we could use some built-in tricks to speed this up.\n", "\n", "First, let's take the PopII parameter set and make a ``ParameterBundle`` object:" ] }, { "cell_type": "code", "execution_count": 4, "id": "b1d2a6a2", "metadata": {}, "outputs": [], "source": [ "popII = ares.util.ParameterBundle(**pars)" ] }, { "cell_type": "markdown", "id": "456624e3", "metadata": {}, "source": [ "This let's us easily extract parameters according to their ID number, and assign new ones" ] }, { "cell_type": "code", "execution_count": 5, "id": "36609007", "metadata": {}, "outputs": [], "source": [ "popIII_uv = popII.pars_by_pop(0, True)\n", "popIII_uv.num = 2\n", "popIII_xr = popII.pars_by_pop(1, True)\n", "popIII_xr.num = 3" ] }, { "cell_type": "markdown", "id": "8491c35d", "metadata": {}, "source": [ "The second argument tells *ARES* to remove the parameter ID numbers.\n", "\n", "Now, we can simply reset the ID numbers and update a few important parameters:" ] }, { "cell_type": "code", "execution_count": 6, "id": "352eeaa0", "metadata": {}, "outputs": [], "source": [ "popIII_uv['pop_Tmin{2}'] = 300\n", "popIII_uv['pop_Tmax{2}'] = 1e4\n", "popIII_uv['pop_temperature{2}'] = 1e5\n", "popIII_uv['pop_fstar{2}'] = 1e-4\n", " \n", "popIII_xr['pop_sfr_model{3}'] = 'link:sfrd:2'\n", "popIII_xr['pop_rad_yield{3}'] = 2.6e39" ] }, { "cell_type": "markdown", "id": "b77af893", "metadata": {}, "source": [ "Now, let's make the final parameter dictionary and run it: " ] }, { "cell_type": "code", "execution_count": 7, "id": "b2e59cf8", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "# Loaded $ARES/input/inits/inits_planck_TTTEEE_lowl_lowE_best.txt.\n", "\n", "############################################################################\n", "## ARES Simulation: Overview ##\n", "############################################################################\n", "## ---------------------------------------------------------------------- ##\n", "## Source Populations ##\n", "## ---------------------------------------------------------------------- ##\n", "## sfrd sed radio O/IR Lya LW LyC Xray RTE ##\n", "## pop #0 : fcoll yes - - x x x - ##\n", "## pop #1 : sfrd->0 yes - - - - - x ##\n", "## pop #2 : fcoll yes - - x x x - ##\n", "## pop #3 : sfrd->2 yes - - - - - x ##\n", "## ---------------------------------------------------------------------- ##\n", "## Physics ##\n", "## ---------------------------------------------------------------------- ##\n", "## cgm_initial_temperature : [10000.0] ##\n", "## clumping_factor : 1 ##\n", "## secondary_ionization : 1 ##\n", "## approx_Salpha : 1 ##\n", "## include_He : False ##\n", "## feedback_LW : False ##\n", "############################################################################\n", "# Loaded $ARES/input/hmf/hmf_ST_planck_TTTEEE_lowl_lowE_best_logM_1400_4-18_z_1201_0-60.hdf5.\n", "# Loaded $ARES/input/hmf/hmf_ST_planck_TTTEEE_lowl_lowE_best_logM_1400_4-18_z_1201_0-60.hdf5.\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ "gs-21cm: 100% |#############################################| Time: 0:00:04 \n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAEsCAYAAADTvkjJAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAAA51ElEQVR4nO3dd3wVZfb48c9JgFACCEovIk0IUqRZUAFdCyKLZVHAxYooi64ulpXVVde2+0Xd9ScquyAiZUVBQUQEpYqgNKWDdJDQSyAEUm/O74+50Syk586dW8779ZpXbmbmzjkMNydPnnnmGVFVjDHGhL8YrxMwxhgTGFbQjTEmQlhBN8aYCGEF3RhjIoQVdGOMiRBW0I0xJkJYQTfGmAhhBd0YYyKEFXRjjIkQVtCNCRARmSwimsdyyOvcTHSwgm5M4DwCtMm1TPGvf8OzjExUEZvLxZjAE5F3gcHAk6r6utf5mOhQxusEjIk0IvIBcBfwB1Ud6XE6JopYQTcmgERkMvA74D5V/cDjdEyUsYJuTICIyAygJ9BfVT/yOh8TfawP3ZgAEJHZwPXAk8BXuTYdUtWD3mRloo0VdGNKSUQEyM5n81RVvS2Y+ZjoZQXdGGMihI1DN8aYCGEF3RhjIoQVdGOMiRBW0ANMRD4TkTT/sltEqopIYxE5KiIZ/q+NPMqtqoikiEiqP7+F/vVByU9EOotIkoik++N/Gqz4IrJFRLJFJC3XuoUi4vOfj1QReS4IMZf7//2pIrJPRM4PQsyQ+PzlJ6+fGa9zCjYROV9EEv2fjXQReaAkx7GCHkAi0gFnHHIdVS2Pc37fACYBy1W1HLAc8GqMcjLQRFUrANWAjiJyfxDzSweGqGoccAHQS0R6BSn+W8Dv81j/papW8C8vBiHmVKCK///gZ+DjIMQMlc/fWQr4mYk284C5/p+N6sDXJTmI3VgUeAJUE5HTQDlgG85t4J38258GVniRmDpDmnLGRFfE+eFR4GKCkJ+qrgHW+F/vF5GjQEIw4qvq2yJyRaCPW9yYqvqPXN8uBO5wOyZB+v8thbx+ZqKGiNQDzgeaAajqKeBUSY5lLfQAUtUfgOnAdiAVOO3/AS7rL2Y5Rc2zX6QiUlZEUoEjwFpVfd+L/PxFpwYw0Yv4ufTwd39s8aAr4h5gVhDihMzn70wF/MxEkytw/nrdKiKnReQnEalRkgNZQQ8gf0HoDrQEKgHl/bPuhQxVzfT/ud8IuFBEbg52DiJSC+duyjdUdW+w4+cyBKgAxAOHKOGfuSUhIl/j3Iw0JFgxQ1E4/MwEQTmcf/urqloR5xfbpyU5kBX0wHoEOKyqP6lqKvA5cBWQKSJtAfxfszzMEQBV3Q2sAu4niPmJSAVgA05/4Z/9qz05P6q6wf8LzofTFdEwGHFFZBRwGdBag3NnX8h9/nLJ72cmmqwGfP6/lgHeA1qU5EBW0ANrA3C+iJzrvx38GmATzn9Yzp+R/8AppEEnIi1yRlWISDWgA/BjsPLzn5P1wB5V7Z1rU1Di55FP21zfPsGv1xfcjPkMcDfQSVWPuh3PbzUh8PnLR34/M1FDVdcBp0TkBv+qPsCukh7MlgAuOBe60oE0YAdQGWgKHAMy/F8be5TbbcBpnD/p0oB5/vVByQ/ngQ/qj5+zPBeM+MBuwOePnwWM9f//pPnzOAC0DULMDP/rnH//hiDEDInPXwE5n/Uz43VOHpyDO3AuhKYC+4FGJTmOzeVijDERwrpcjDEmQlhBN8aYCGEF3RhjIoQVdGOMiRBhWdBFZJeIrBOR1SKy0r+uuojMEZGt/q/VvM4TQEQmeJ1DYbzMMdpiR0vM4gqHHN0WiHMQlqNcRGQX0FFVj+RaNxw4pqr/EJGngWr6640rnhGRU6paqYTvfVFVC50BsKj7FfD+EudYWtEWO1piFlc45Oi2QJyDsGyh56M3MM7/ehxws3epBEyfAO9njIlg4dpC3wkk4dw88R9VHSUix1X1nFz7JKnqWd0uIjIIGOT/tkNMjLu/07KzsylpjKK+tzQxAvH+0oi22NESs7jCIUe3FefnXVUlr23hWtDrquo+EakJzMGZD+LzohT03GJjY9Xn87mbrDHGBJCIZKtqbF7bwvJXoqru8389BEwDOgMHRaQOgP/rIe8yNMaY4Au7gi4ilUSkcs5r4DqcCZ8+x5n0CP/X6d5kaIwx3giZie6LoRYwzZmYjTLAh6o6W0RWAJP9j1T7GbtQaIyJMmHZhx4o1odujAk3EdeHbowx5mxW0I0xJkKEYx96qYlIL6CXvx/eGGMigvWhWx+6MSaMWB+6McZEASvoxhgTIaygG2NMhLCCbowxEcIKujHGRAgbtmiMMRHChi3asEVjTBixYYvGGBMFrKAbY0yEsIJujDERwgq6McZECCvoxhgTIWzYojHGRAgbtmjDFo0xYcSGLRpjTBSwgm6MMRHCCroxxkQIK+jGGBMhrKAbY0yEsGGLxhgTIWzYog1bNOZ/+Hxw/DicOAGpqZCW5nzNWdLSnMXng+zs/BdViI39dSlTpuDX5cpBXFzhS7lyEM1tsYKGLUZlC92YaJOVBbt3Q2Ii7NsHe/c6X/ftg4MHISnJKeJJSZCc7HW2hSuo4Jcv7yw5r8/8GqhtcXEQE2Kd1lbQjYkgKSmwdi2sXg2bN8PWrc6ya5dT1HOrWBHq1oVataBBA2jTBqpV+3WpWhUqVDh7ySluZco4BS2/BZyWelaW05r3+c5+nZWlnDyZyvHjJzlxIpWTJzM4eTKDlJRMUlIyOXUqi9Onszh9Opu0NCUzU8jIiCEjQ8jMjCEzU8jMjCUrK4a0tFhSUmLx+cr4l3L+r2X/52t2dlmyssoQiEuIsbE+ypb1UbZsNmXL+ihXLpuyZbMpV04pVy6buDj1F3894xeDUKFCzhJDhQpCxYoxVKgQQ8WKsVSqFEOlSmUoX17O+qVSkIjqchGRG4D/B8QC76nqPwra37pcTDjz+ZzivWgRfP+9U8S3bHG6OgAqVYJmzX5dmjaFhg2dIl63LlSp4k7XRWZmJvv372ffvn3s3buXffv2/bIcOXKEpKQkjh07RlJSEklJSWRmZpYoTtmyZSlTpsz/LLGxscT4f5vk1La8vqqCahlUywHl8fnKkJV1dvGHOP9S3r/EnfE1r3XF3b+48u9yiZiCLiKxwBbgWiARWAH0U9WN+b3HCroJN9u2wRdfwJw5sHjxr90jDRpA+/Zw8cW/LvXru9fXrKokJiaydu1aNm/ezLZt29i+fTvbtm1j9+7dnPlzVaZMGerUqUONGjWoXr061apV+2XJ+b5y5cpUrFjxl6VSpUq/vK5QoQJxcXH/U7iDITs7G5/Ph8/nIysrq8hfcy+ZmZlnfc15nZGRSVpaNqmpSmqqkpbmfE1Pd65TpKcL6emQkSH+18L8+fdHRUG/DHhBVa/3fz8MQFX/nt97YmJidPnyvTRufA7Vq1cIUqbGFJ0qrFwJkyc7hfynn5z1F14IXbvCVVc5S4MGbuag7N69myVLlrBixQrWrFnDmjVrSEpK+mWfc845h6ZNm9K0aVOaNGlCo0aNqFevHnXr1qVu3bqce+65v7ScTelEy0XResCeXN8nApcU9AZVoVOnOoCPyy9fyLffdiUmJoovn5uQsXMnTJzoLFu2QNmy0L07/OEP0LMnNG7sbvzt27cze/ZsvvnmG5YsWcK+ffsAqFixIq1bt6ZPnz60adOGtm3bkpCQQPXq1d1NyBRJJBX0vCrxWX9+iMggYFDOW/r3X8S335bhu++6ccstC5k+vZubORqTr+xs+PprGDECZs1yWufdusFTT8Ftt8E557gXOysri4ULF/LFF18wa9YstmzZAkCDBg3o2rUrXbp04YorruCiiy4KWneHKb6o7nLJ6UPPysqmYcPl7N/fnqlTd3PLLc2ClLUxkJkJ48fD8OFOa7xWLXjwQbj/fucipluys7P57rvvmDRpElOmTOHw4cOUL1+ebt26ceONN9KjRw+aNm3qXgKmRArqcomkgl4G56LoNcBenIui/VV1Q37vyX1RdNOmI7RqpVSuvI+kpDbW9WJcl1PIX37ZGVbYoQP86U/Qp49z84xbjh49ytixY/n3v//N9u3bKV++PL169aJfv35cf/31VKxY0b3gptQKKuj+ITyRsQA34hT17cAzhe0fExOjufXt+42C6vPPL1dj3JKdrTpjhmqzZs7guU6dVL/4wlnvpo0bN+o999yjcXFxCugVV1yh48eP1+TkZHcDm4ACfJpPTYuYFnpJnDls8fTpTKpW3UeZMmmcOtU8IlrpR44cYf78+SxduoHDhw9SrVocl156Cb1796ZSpUpepxd1Nm1yWuFffQUtWsBrrzkXOd28lX3VqlW8+uqrfPrpp1SoUIG7776bwYMH07p1a/eCGtdETQu9uMuZLXRV1YEDv1VQ/dvfwreVnpWVpWPGzNYWLd5T+FrhuDqX2FRFDir8R887r7POmDHD61SjRmam6iuvqJYrp1q1quqbb6pmZLgbc8eOHXr77bcroFWqVNFnnnlGDx065G5Q4zqshZ63vG4sSknJ4JxzjhAff4jjx9t5k1gJZWdn8/LL83ntNSElpRsQS40aB7jySujcuQYisaxapUyblk1GRgaq9zJxYi/uvPNOr1OPaOvWwb33wg8/wO23O6NYatZ0L15ycjKvvPIKb775JrGxsTzxxBMMHTqUc9wcJmOCxlroZyxAL2CUiOT5G/C3v12goPree+uK/FvTa//613caH79YQTU29qjefPNG/emnzDz3TUxUvfzyLIUsjY39nf74449BzjY6ZGervvOO0yqvUUN1yhT3Y37++edar149FRG9++67NTEx0f2gJqgooIXueXH1csmry0VVde/eZBVJ0rp1vy/6WfbIypV7tGHDef5CfkT79/9RT5zIKvR9KSmq7dtnqkiyJiT00Kyswt9jiu7ECdU77nB+wnr0UHW7p+Pw4cPat29fBbR169a6bNkydwMaz1hBL2ZBV1W9/PIFCj79+uudRTrJXnjwwUUKSQrp2rXrEj1yJK1Y79+xQ7V8+QyFGfrhhx+6lGX02bpV9cILVWNiVF99VdXnczfeggULtG7dulquXDl98cUXNT093d2AxlNW0EtQ0NesOaiQpgkJ3xTtLAfRgQMntWlTZ4hllSprdN68kv9Z/eqrPgXV88//vfrcrjxRYNEi1erVnWXBAndjZWVl6QsvvKAxMTHavHlzXbVqlbsBTUiwgl6Cgq6qeuGF3yik6vr1oTMy4JNPftJy5bYq+PTKKxdqWlrpukpOnVKtWvW0wgydM2dOgLKMTv/9r9Nf3ry500p304kTJ7RHjx4K6IABA/TkyZPuBjQhwwp6CQv6zJnbFVSvvHJBoSc5GO6/f5HCaY2JOajDhwfuQubTT2cq+PSmmx4J2DGjzciRzk9T166qR4+6G2vHjh2akJCgZcqU0X//+9/uBjMhxwp6CQu6qmrt2ktV5KgePJhS6L5uSU/P0vbtFyioVqv2o65dG9i/GPbsURXJ0tjY1/TEiRMBPXY0+Oc/nZ+knj1VU1PdjbVy5Uo977zztFq1ajpv3jx3g5mQVFBBj8oJikWkl4iMcs5NwZ57rgKq1RkyZGUQMjvbgQMp1K//Az/+2I3WrRexd29rWreuEdAY9evDZZedwOe7nS+/nB3QY0e6v/8dhg6F3/0Opk4t/BFhpbFo0SK6d+9OpUqVWLp0KVdffbV7wUx4yq/SR8NSlBa6qmp8/FotU2a3pqbmPa7bLTt3Htf4+HUKWdq3r7sXZ8eOdS6OXnvtX12NE0lGjHBa5v37O3eCumnWrFlavnx5bdGihe7Zs8fdYCakYS300vnDH06RldWQYcOWBy3mhg2HadnyACkpzXniieVMmnSVq/FuvjmGmJgsvvmmJhkZGa7GigQffgiPPAK//S2MG+c8MNktc+fOpXfv3rRs2ZJFixZRv35994KZsGa3/hfhmaIZGT7i43+mbNlUTp5s6fqkXStW7OeKK1LJyKjNyy9v5JlnOroaL0f79gdYteo4ixcfpUuXLkGJGY6+/BJ694Yrr3Reu9nN8t1333HttdfSuHFjvvnmG3sykCnw1n9roRdBuXKx3HlnIqdPJzBs2DJXY82fv5vLLvORkXEub7+9NWjFHOCWWyoDLfjssx+CFjPcrF7tzMfSti189pm7xXz16tXceOON1K1blzlz5lgxN4WyFnoRWugAaWlZVK26G4ATJ86nfPnA/409ffo2br01HtWyjBt3kAEDEgIeoyDr1kGbNtCy5ets3PhEUGOHgwMHoHNnZ97KFSugdm33Yu3du5fOnTsTExPD4sWLOf/8890LZsKKtdADoHz5Mjz22GEyMprw4IPfB/z4Eydu4pZbqgEwbVpS0Is5wEUXQcWKJ9iypQHp6elBjx/K0tLg5pvh6FH4/HN3i/np06fp3bs3ycnJzJw504q5Kbr8rpZG8kIhsy3mx+fL1vj4tRoTs1937z5erPcWZMSI1QrHNTb2Z503b1fAjlsSXbv+rHBQFy361tM8Qs299zojWj791N04Pp9Pb7vtNhURm6/e5Akb5fK/VHWGqg6SYj4mJiZGeOedWLKza3D99WsCkstzzy3nkUeaUa7cUb7/vgxXX+1ta6x376pATWbM2OxpHqHkgw9g7Fh49lm49VZ3Y+U8Wej111/npptucjeYiTz5VfpoWIo6Dv1MHTsuUFB9/fXS3X7fv/83CllaseKGkJkvZsMG9T/n8m2vUwkJ69apVqig2r27qtszDM+fP19jYmK0f//+mu32A0ZN2MKeWJS34lwUze3IkdPUq3eArKx4VqzIpn374nWoZmcrV131DUuWdKNGjRWsX59AzZqh8XzP7GwoXz6FuLjPOXmyv9fpeColBTp2hOPHndEtbvabHzx4kHbt2lG1alVWrlxJfHy8e8FMWLOLogF23nkVmTw5m+zsinTtepgjR04X+b0//3yC889fypIl3bjwwkX8/PPFIVPMAWJioEmTw6SktObgwYNep+Opxx+HLVtg0iR3i7nP5+POO+/kxIkTTJkyxYq5KTEr6CXUu3dThg5dR0pKK5o128yBAymFvmfs2A00aXKCxMRO9OixkI0br3Rl+GNpdekSA7Ri/vzoHY8+cyaMGgVPPgndu7sb680332TevHmMGDGC1q1buxvMRLb8+mKiYSlpH3pugwcvVsjSMmV26dNPL9XMzP99SITPl60ff/yTtm69UCFTY2P36H/+s7bUcd00Y8ZpBdUBA8Z7nYonDh9WrVVLtXVr1bTiPQSq2DZu3KhxcXHau3dv6zc3RYL1oeetpH3oZxoxYg2PP16FzMwLEDlMrVo7qFgxk+TkOJKS6uHz1QV8tGq1hBkz2nDBBeeUOqabjh2Dc8+FVq0msH79AK/TCbo+fWD6dOfmobZt3YuTlZXF5Zdfzo4dO9iwYQO1atVyL5iJGNaHfobiTJ9bFI880pbk5AY89th3nH/+FpKSqrN7dyNOnYqnbt1d3H33t6xZc5T1668K+WIOUL06VKx4iN27q3mdStBNmwaffAIvvuhuMQcYPnw4K1as4N1337VibgLCWugBaKFHooSEn9i0SUhOrkvlypW9TicoTpyAhASoUcNpnZct616szZs306ZNG3r37s3kyZPdC2QijrXQTbG1awfQjGXLNnicSfA88wzs3w+jR7tbzFWVhx9+mAoVKvDWW2+5F8hEHSvoJk+/+U11IIavvjrgdSpBsXQpvPuuM8d5p07uxpoyZQpz587llVdeobab4yFN1AmrLhcReQF4ADjsX/UXVf3Sv20YcD/gA/6oql8VdjzrcslfYqLSoIHQpcvHLF58h9fpuCozEzp0gKQk2LgR3OxhSk5OpkWLFtSpU4fly5cTG5vnX87G5KugLpfQGwRduH+p6uu5V4hIAtAXaAXUBeaKSHNVtWpdQvXqCWXLJrF1a+jc9OSWkSOdqYOnTnW3mAO88MILHDhwgM8++8yKuQm4SOly6Q18pKrpqroT2AZ09jinsCYCtWsf4MiRBoTTX3HFdfQovPAC/OY3zvS4btq6dSsjRoxg4MCBdO5sH08TeOFY0B8WkbUi8r6I5IyrqwfsybVPon/dWURkkIisFJGVkVyoAqFFi1Sys1uya9der1NxzfPPQ3Iy/Otfzi8xNw0bNoy4uDheeukldwOZqBVyBV1E5orI+jyW3sBIoAnQDtgPvJHztjwOlWe1VtVRqtpRVTsWd/rcaNOxY3mgHPPm/ex1Kq5Yv97pbnnoIefhHm767rvv+PTTT3nqqadszLlxTcgVdFX9japelMcyXVUPqqpPVbOB0fzarZIINMh1mPrAvmDnHmmuvromAEuWnPA4k8BThT/9CapWhb/9ze1YyhNPPEGdOnV4/PHH3Q1molrIFfSCiEidXN/eAqz3v/4c6CsicSJyAdAMWB7s/CJNly7nAj7Wrcv2OpWA+/JLmDvXKebnnuturKlTp/L999/z4osvUqlS5F9kNt4Jt2GLE3C6WxTYBTyoqvv9254B7gOygMdUdVZhx7Nhi4WrUGE38fG7OXz4Kq9TCZjsbOfGqdRUZ5iimzcR+Xw+LrroImJjY1mzZo2NbDGlFjHDFlU135miVPUV4JUgphMVatU6TGJiXa/TCKhJk5xhipMmuVvMAT7++GN++uknPvnkEyvmxnXFbqGLyIVAY6ACzg0+q1S18MnAQ5C10At37bXfMXfuJSQmnqBevepep1NqGRnQooXTd/7DD84DPdzi8/lo1aoV5cqVY/Xq1cS4GcxEjVK30EWkEfAH4PdALf53VEmWiCwG/g1M0TDowxGRXkAvG+VSuA4d4pg7N5avv/6Ze+8N/4L+3nuwc6fTh+52ff3444/ZvHkzn3zyiRVzExSFttBF5DVgMDAPmAYswxlVkgZUB1oD3YD+wEngXlX90b2UA8da6IX76qufueGGhtx330LGjOnmdTqlcuoUNGkCF14ICxe6O+7cWufGLaVtoVcGmqtqXsMAD/qXucCzItIHaAmERUE3hevatR6QHhEjXf7f/4ODB51b/N3+48xa58YLYTXKJdCshV40FSpspnLlExw6FL63qycnw/nnwxVXwIwZ7sZSVdq0aQPAmjVrrKCbgCr1fOgiEleEfdoVMy8TJmrWPEJSUniPdHnnHTh+3Jm3xW2zZ89m/fr1PPXUU1bMTVAV9dP2YUEb/bMdFjpdrQlPTZumkZVVn4MHT3udSomkpMAbb8CNNzrT5Lpt+PDh1K9fn759+7ofzJhcilrQLxeREXltEJGmOH3odmdmhGrfvhwAc+aE5yRdI0c6syr+9a/ux1q+fDkLFy7kT3/6E2XdHuRuzBmKWtB7Anf7HyLxCxFpiDP6ZRNwW4Bzc02gHxId6bp3rwHAt98meZxJ8Z0+Da+/DtddB5de6n681157japVq/LAAw+4H8yYMxSpoPuHIf4OeF5EBsAv86rMA/YCvVQ1w7UsA0xVZ6jqIBuHXjRdu54PnGLNmvAb6TJqFBw6FJzW+bZt25g6dSqDBw+Omgdrm9BSrFEu/mI+GudRb3/BGYt+taqG5XR8Nsql6OLi1lCtmnDgQBuvUymytDRo3Ni5M3T+fPfjPfzww4wePZpdu3ZRp06dwt9gTAkEbC4XVZ0gIrWB8TjdLNeGazE3xVOz5kEOHGjvdRrFMmYM7N8PHxZ4ST8wkpOTGTduHHfccYcVc+OZot76//UZqzL9y0e5uy1U9brApWZCSZMmqSQmnseBAz5q1w79SaaysuC116BLF+ja1f1448aNIyUlhUceecT9YMbko6gXRfeesUwCVuex3kSoiy92RmwsWHDI40yKZsoU2L0b/vxn9+8Kzc7O5u233+aSSy6hU6dO7gYzpgBFaqGr6r1uJ2JCW9eu5/Lmm7BoURL9+oV2l4IqDB/u9J337Ol+vLlz57JlyxYmTJjgfjBjChBW86Eb73Tp0gQ4yurVWV6nUqh582D1aqcPPRg3ao4YMYKaNWvSp08f94MZU4BiF3QRuR64BqjJGV02qnpXgPJylU2fW3w1apxHmTJL2Lkz9B9wPHw41KkDd97pfqwdO3Ywc+ZMnn32WeLiCp0hwxhXFav9IiIvA7OA64DaQI0zlrBg49BLpkaNAxw5UptQvh9r9WqYMwcefRSCUV/fffddYmNjeeihh9wPZkwhittCHwTco6rj3UjGhLbGjVPZvz+ePXugYUOvs8nba69BfDw8+KD7sdLS0hg7diy33HILdeuG9+RlJjIUt4cxG/jOjURM6GvXzhmuuHhxaN56sHs3fPyxU8zPOcf9eNOmTePYsWMMGjTI/WDGFEFxC/q7wEA3EjGh76qrzgVCd06Xf/3LGaL46KPBiffee+9xwQUXcPXVVwcnoDGFKG6Xy0vAFyKyBliLc3PRL1T1vkAlZkJPx45NgT0hOdLl2DHneaH9+0ODBu7H2759O/Pnz+fll1+2Oc9NyChuQX8R6AFsBOoAIXx5zATa+eefj8hctm+/yOtUzjJypPPM0CeeCE689957j9jYWO69127RMKGjuJNzJQF/UtUPXMsoCHINW3wgOzv8ZhD0Uo0aYzl27Pekp5elTIjcxZCW5jxerkMH+PJL9+NlZmbSsGFDOnfuzPTp090PaEwupX4EXS4ZwOLSp+QtG7ZYco0bnyY7uyzbtnmdya/Gj3emyH3yyeDEmzlzJgcOHGDgQLucZEJLcQv6KJypc02UatvWaRj88ENoTH/v8zmPl+vYEbp1C07M0aNHU7duXXr06BGcgMYUUXELeh3gIRH5UUTGisio3IsbCZrQ0qVLdSCLb7457nUqAHz+OWzZ4rTOg/EH1549e5g9ezb33XcfZUKlz8kYv+IW9CY4syyeABoBzXItTQOVlIj0EZENIpItIh3P2DZMRLaJyGb/NAQ56zuIyDr/trfE+lNc0bZtc2ATK1d6/2CQnEm4LrgAbr01ODHHjh2LqnL//faHqgk9xX3ARXe3EjnDeuBW4D+5V4pIAtAXaAXUBeaKSHNV9QEjce5kXQp8CdyAM02BCaDmzZsDU9i6NUgVtABLlsDSpfD22wTlAq2q8sEHH3DNNdfQqFEj9wMaU0yFttBFpENRDyYi5UWkZelSAlXdpKqb89jUG/hIVdNVdSewDejsf75pFVX9Xp1hO+OBm0ubhzlbxYoVqVZtFykpldm/39tcXnsNzj0XgjVycPHixezcuZO77747OAGNKaaidLlMF5FpInK9iOS5v4jUE5FhwFagS0Az/F/1gD25vk/0r6vnf33m+rOIyCARWSkiK4szZNP8qlmzkwD8+KN3Ofz0k9N/PmQIVKwYnJjjx48nPj6eW265JTgBjSmmovyheiHwNDARKC8iq3CeTpQGVMfp/rgAWAj0U9UiDWsUkbk4Mzae6RlVzW9wb1794lrA+rNXqo7CGa1DbGysVfQS6NixDMuXww8/ZNOzpzd3Sb7xBpQvDw8/HJx4qampTJ48md/97ndUqlQpOEGNKaZCC7qqngL+6p86twdwFU4BrwYcBt4BZqnqT8UJrKq/KX66JAK5b+yuD+zzr6+fx3rjgjZtLgA2s2RJAyBIzeNcDhxwxp7fdx/UCNKkzdOnTyc5OZm77gqLKf9NlCrypSRVTQc+8y9e+Rz4UET+iXNRtBmwXFV9InJSRC4FlgF3ASM8zDOitWzZEviRVau8mUN3xAjIzIShQ4MXc/z48TRo0ICuwXjitDElFJKzConILSKSCFwGzBSRrwBUdQMwGWcumdnAEP8IF4DBwHs4F0q3YyNcXNOiRQtgBYcPV+DAgeDGTkmBd9+FW26BZs2CE3P//v189dVXDBgwwCbiMiEtJO+MUNVpwLR8tr0CvJLH+pVA6M0aFYFq1KhB5cobOHkSvv/eKa7B8t57cPw4PPVU8GJ++OGHZGdnW3eLCXnW3DDFJiK0apWOSAbfBfFxJ5mZzpznV14Jl1wSnJiqyrhx47jkkku48MILgxPUmBKKyoIuIr1EZJQNWyy5hIQmlCmzlu+/D17MyZPh55+DNwkXwJo1a1i3bp21zk1YiMqCbrMtll6LFi3IzPyGlSuVjCDM05WdDa++Cq1aQc+e7sfLMX78eMqWLUvfvn2DF9SYEorKgm5Kzxnp8j3p6cLq1e7HmzYNNm6EZ5+FYF2XzMzM5L///S+9evWievXqwQlqTCmU6kdDRK7zz7y4W0Smi0j7QCVmQltCQgLg9LcsWeJuLFV4+WVnVEufPu7Gyu3rr7/m0KFD1t1iwkZp2zojgUdxRpf8C/iniAwodVYm5DVq1IgqVVKoWvUQ8+e7G2vmTFi9Gv7yF4jN8zkt7hg/fjznnXeezXtuwkZpC/ohVf1WVU+q6kLgRuDPpU/LhLqYmBjatm1L+fKLWbjQGYHihpzWeaNGcOed7sTIS1JSEtOnT6dfv36UK1cueIGNKYUSFXQR+VhEngKWiMhLIlLWv8kHpAcsOxPS2rVrx/Hjn5CSAsuWuRNj7lzn2E8/DWXLFr5/oEyZMoX09HTrbjFhpaQt9BHAaeAcnPldtovIQmALzlzkJgq0a9eO9PQviYlR5swJ/PFV4a9/hfr14Z57An/8gowbN46EhAQ6dCjy7NHGeK5IBV1ERovIL7MwqepiVX1bVQeqakecyboeBp4FyruTauDYOPTAaNu2LXCCJk2OuVLQp093WucvvABxcYE/fn62bdvGd999x4ABA7ChrSacSFGKmoj4gDqqesj9lIInNjZWfT7vH6UWrtLS0oiPj+fSS+eydGk3Dh+GatUCc2yfD9q0cb6uXx+cJxLleP7553nppZf4+eefqV+/fuFvMCaIRCRbVfMcHlDULhdrppizlC9f3j8e/XN8PpgxI3DHnjjRGXf+8svBLeaqyoQJE7jmmmusmJuwU5w+dOufMGdp164du3ZNoUED+PTTwBwzNRWeew46dIDbbgvMMYtqyZIl7Ny50y6GmrBUnIL+HxF5UkS6i0gV1zIyYaVdu3bs3ZtIjx6n+eorOHmy9Mf8v/9z5mx5/XUIdhf2+PHjqVixoj1mzoSl4hT0ujgXPecBSSKyRUQ+FJGhItJVRCq7k6IJZR07dgSgSZPVpKeXvpW+c6dT0Pv2hW7dSp9fceQ8Zu62224jPj4+uMGNCYDiFPTf4gxTTMB5ItBMnMfB/Q1YACQFOjkT+jp06EBMTAwnT86meXMYM6Z0xxs61Lkb9LXXApNfccyYMYMTJ05Yd4sJW0W93KQA6gyJ+cm//BdAnHFdCUDYzOMiIr2AXjYkrfTi4+O56KKLWLFiOQMHOg+e2LQJWrYs/rGmToXPPoO//90Zex5sEyZMoF69enTv3j34wY0JgFKPclHHBlWdEKCcXGfT5wbWJZdcwvLlyxkwQClXDt58s/jHOHQIHnzQuRD6+OMBT7EI8Q8xa9Ysfv/73xMbzAljjAmgohb0a4DjLuZhwtgll1xCUlISyclbue8+GDsWEhOL/n5VGDTIuaA6blxwb/HPMWnSJHw+HwMG2NxyJnwVtaA/DPQGEJEWIvKViKwXkY9E5Eb30jPh4BL/8+CWLVvGn//sFOgXXyz6+4cPd+4K/fvfnQdYeGH8+PG0b9+eVl4lYEwAFLWgdwXW+l9/CBwDvgUaA1+IyPti/RdRq2XLlsTHx7Ns2TIaNYI//hFGj6ZIj6ebPt2ZFveOO+Cxx9zONG8bNmzgxx9/tIuhJuwVtaDHA8ki0gZ4T1X7qepgVe2MU+x7AIPcStKEttjYWDp16sTSpUsB+NvfoEEDGDAAjh7N/32zZ8Pttzv95mPGBH/MeY4JEyYQGxtLv379vEnAmAApakE/CNQBugOTc29Q1W+Bx4AHA5qZCSuXX345q1ev5uTJk8THOw90TkyEm25yLnjmlp0NI0ZAr16QkABffQWVKnmTt8/nY+LEidxwww3UrFnTmySMCZCiFvQvgdHAk0BenYw/AE0DlZTbbLbFwOvWrRs+n48l/ufRXXopTJrkPGmobVv4xz9g1iwYORI6d3a6Za6/HhYuDNyEXiUxZ84c9u7dy9133+1dEsYESFFnW6wCvIEzNe5inGGME1U1xb/9eaC/ql7oYq4BZ7MtBs6pU6eoVq0aQ4cO5R//+Mcv61etcm4WWrjw132bN4dnnnG6ZLy+8tKnTx8WLFjA3r17iQvmHL3GlFBBsy0W6cYiVU0GHsh1wL8DR0RkN1AJpzvmTwHI1YSpSpUq0blzZxbmrtzAxRfDggWwd69zW3/t2tCkifeFHODw4cNMnz6dIUOGWDE3EaFETyxS1WFAO2AcTp/67ar6VgDzMmGoW7durFy5kpN5zNBVrx5ccQU0bRoaxRxg4sSJZGZmcv/993udijEBUeKHRKvqT6r6qqoOVdUATZzqEJE+IrJBRLJFpGOu9Y1EJFVEVvuXf+fa1kFE1onINhF5y4ZRBt+Z/eihTFUZM2YMnTt35qKLLvI6HWMCosQF3WXrgVuBRXls266q7fzLQ7nWj8QZOtnMv9zgfpomt8suu4yyZcsyf/58r1Mp1PLly9mwYYO1zk1ECcmCrqqbVHVzUfcXkTpAFVX93j+B2HjgZrfyM3mrVKkSV155JV9+GfrPCR8zZgwVKlSgb9++XqdiTMCEZEEvxAUiskpEvhGRK/3r6gG5Zw9J9K87i4gMEpGVIrLShi0GXs+ePdmwYQO7d+/2OpV8nTp1io8++og+ffpQpYo9q8VEDs8KuojM9c8Hc+bSu4C37QcaqurFwFDgQ/+Qyrz6y/Os1qo6SlU7qmpH62YPvJ49ewKEdCt90qRJnDx5koEDB3qdijEBVaRx6F4RkYXAE6q6sqDtwF5ggaq28K/vB3RT1QLvXrVx6IGnqjRt2pSWLVvyxRdfeJ3OWVSV9u3bk5WVxdq1a7Ff6ibcFDQOPay6XESkhojE+l83xrn4uUNV9wMnReRS/+iWu4DpHqYatUSEnj17Mn/+fFJTU71O5yzLli1j9erVDBkyxIq5iTghWdBF5BYRSQQuA2aKyFf+TVcBa0VkDfAJ8JCqHvNvGwy8B2wDtgOzgpy28bvppptITU1lzpw5XqdylnfeeYfKlStz5513ep2KMQEX0l0ubrMuF3dkZGRQu3ZtbrzxRiZOnOh1Or84fPgw9evXZ9CgQYwYMcLrdIwpkYjpcjHhoVy5ctx6661Mnz49pLpdxowZQ0ZGBoMHD/Y6FWNcYQXduOKOO+4gJSWF2bNne50KAJmZmbz77rt0796dhIQEr9MxxhVRWdBt+lz3de/enfPOO4+PPvrI61QAmDx5Mnv27OFxL55AbUyQWB+69aG75pFHHmH06NHs27eP6tWre5aHqnLxxReTmZnJunXriImJynaMiRDWh248MXDgQNLT05kwYYKnecybN481a9bw+OOPWzE3Ec1a6NZCd1Xnzp1JTU319Cae66+/nrVr17Jr1y6b99yEPWuhG8888MADrF+//pcHSAfbihUr+Prrr/njH/9oxdxEPGuhWwvdVSdPnqRBgwZcd911TJ48ufA3BFjPnj1ZunQpO3futIm4TESwFrrxTOXKlRk8eDCffvop27dvD2rsZcuW8eWXX/LEE09YMTdRISpb6CLSC+glIg9kZ2d7nU7E279/P40aNWLgwIG88847QYvbo0cPVqxYwc6dO6lcuXLQ4hrjJmuhn0FVZ6jqIJucKTjq1KnDgAEDeP/999m3b19QYs6fP5/Zs2fz1FNPWTE3USMqW+g5rA89eHbs2EGLFi245557GDVqlKuxfD4f7du3Jzk5mU2bNlG+fHlX4xkTTNZCN55r3LgxgwcPZsyYMWzatMnVWO+//z5r165l+PDhVsxNVLEWurXQg+bw4cM0adKErl27MmPGDFdiHD16lISEBJo3b86iRYtsznMTcayFbkJCjRo1+Otf/8oXX3zB1KlTXYkxdOhQjh07xttvv23F3EQda6FbCz2osrKy6NSpEwcOHGDjxo1Uq1YtYMeeNWsWN954I88++ywvvfRSwI5rTCgpqIUelQXdhi1664cffuDSSy/lt7/9LZ988klAWtL79++nXbt21KhRgx9++MHuCjURy7pczmDDFr3VoUMHXn31VaZOnRqQJwdlZWXRr18/UlJSmDx5shVzE7WisqAb7z3++OP06tWLoUOHMnPmzBIfR1UZMmQI33zzDSNHjrSHV5ioZgXdeCImJoYPP/yQtm3bcvvtt7Nw4cJiH0NVee655xg1ahTDhg3jrrvuCnyixoQRK+jGM/Hx8cycOZNGjRpxww03FGvki8/nY+jQobz88svcd999vPLKKy5makx4sIJuPFW7dm0WLVpEu3btuO2223j00Uc5efJkge/Zs2cP1157LW+++SaPPfYYo0ePtiGKxhClo1xy2LDF0JGens6TTz7JiBEjqFWrFo899hj9+vWjYcOGiAiqyrp16/jggw8YOXIksbGxvPXWW9x7771WzE1UsWGLZ7Bhi6Fr2bJlDBs2jAULFgBw3nnnUaVKFQ4dOkRKSgqxsbH079+fF198kUaNGnmbrDEesIKeD2uhh66tW7cye/Zs1q1bx+nTp6levTpt27alV69e1KxZ0+v0jPGMFfR8WEE3xoQbu7HIGGOiQEgWdBF5TUR+EpG1IjJNRM7JtW2YiGwTkc0icn2u9R1EZJ1/21tiV8qMMVEmJAs6MAe4SFXbAFuAYQAikgD0BVoBNwDvikjOnx4jgUFAM/9yQ7CTNsYYL4VkQVfVr1U1y//tUqC+/3Vv4CNVTVfVncA2oLOI1AGqqOr36lwUGA/cHOy8jTHGSyFZ0M9wHzDL/7oesCfXtkT/unr+12euN8aYqFHGq8AiMheoncemZ1R1un+fZ4As4L85b8tjfy1gfV5xB+F0zdgNKcaYiOJZQVfV3xS0XUTuBm4CrtFfx1YmAg1y7VYf2OdfXz+P9XnFHQWMAmfYYomSN8aYEBSSXS4icgPwZ+C3qno616bPgb4iEiciF+Bc/FyuqvuBkyJyqX90y13A9KAnbowxHvKshV6It4E4YI6/W2Spqj6kqhtEZDKwEacrZoiq5twZNBj4AKiA0+c+66yjGmNMBLM7Re1OUWNMGLE7RY0xJgpYQTfGmAgRqn3orso1fa7XqRhjTMBYH7r1oRtjwoj1oRtjTBSwgm6MMRHCCroxxkQIK+jGGBMhrKAbY0yEsGGLxhgTIWzYog1bNMaEERu2aIwxUcAKujHGRAgr6MYYEyGsoBtjTISwgm6MMRHChi0aY0yEsGGLNmzRGBNGbNiiMcZEASvoxhgTIaygG2NMhLCCbowxEcIKujHGRAgr6MYYEyFsHLoxxkQIG4du49CNMWHExqEbY0wUsIJujDERIiQLuoi8JiI/ichaEZkmIuf41zcSkVQRWe1f/p3rPR1EZJ2IbBORt8Q6yI0xUSYkCzowB7hIVdsAW4BhubZtV9V2/uWhXOtHAoOAZv7lhqBla4wxISAkC7qqfq2qWf5vlwL1C9pfROoAVVT1e3Wu8o4HbnY3S2OMCS3hMGzxPuDjXN9fICKrgGTgWVX9FqgHJObaJ9G/7iwiMginJZ/zfXYJ8xKgJEOEgv0+L2KGU65exAynXL2IGU65ehEz/4a4qnqyAHOB9XksvXPt8wwwjV+HV8YB5/pfdwD2AFWATsDcXO+7EphRhBxWliL/UeHwvlLGLNH58ShXOz+hF9POjwvnp6DFsxa6qv6moO0icjdwE3CN+v/1qpoOpPtf/yAi24HmOC3y3N0y9YF9buSdy4wweV9p3xvseHZ+3HmvnZ/QixlwIXljkYjcAPwT6Kqqh3OtrwEcU1WfiDQGvgVaq+oxEVkBPAIsA74ERqjql4XEWamqHV37h4Q5Oz8Fs/NTMDs/BXPj/IRqH/rbON0rc/yjD5eqM6LlKuBFEckCfMBDqnrM/57BwAdABWCWfynMqADnHWns/BTMzk/B7PwULODnJyRb6MYYY4ovJIctGmOMKT4r6MYYEyGitqCLyA0istk/VcDTXucTCkRkl3/6hNUistK/rrqIzBGRrf6v1bzOM1hE5H0ROSQi63Oty/d8iMgw/+dps4hc703WwZPP+XlBRPbmmp7jxlzboub8iEgDEVkgIptEZIOIPOpf7+7nJ9DjIMNhAWKB7UBjoBywBkjwOi+vF2AXcN4Z64YDT/tfPw38n9d5BvF8XAW0B9YXdj6ABP/nKA64wP/5ivX63+DB+XkBeCKPfaPq/AB1gPb+15VxpjBJcPvzE60t9M7ANlXdoaoZwEdAb49zClW9gXH+1+OIoikVVHURcOyM1fmdj97AR6qarqo7gW04n7OIlc/5yU9UnR9V3a+qP/pfnwQ24dy97urnJ1oLej2cu0xz5DtVQJRR4GsR+cE/RQJALVXdD86HFKjpWXahIb/zYZ+pXz3snyn1/VxdClF7fkSkEXAxzj0yrn5+orWg5zW1ro3fhC6q2h7oAQwRkau8TiiM2GfKMRJoArQD9gNv+NdH5fkRkXjgU+AxVU0uaNc81hX7/ERrQU8EGuT6PhhTBYQ8Vd3n/3oIZw6dzsBB/2yWObNaHvIuw5CQ3/mwzxSgqgdV1aeq2cBofu02iLrzIyJlcYr5f1V1qn+1q5+faC3oK4BmInKBiJQD+gKfe5yTp0SkkohUznkNXIczWdrnwN3+3e4GpnuTYcjI73x8DvQVkTgRuQBnTv7lHuTnqZxi5XcLzmcIouz8+B+wMwbYpKr/zLXJ1c9PqN767ypVzRKRh4GvcEa8vK+qGzxOy2u1gGn+qRbKAB+q6mz/HDmTReR+4Gegj4c5BpWITAK6AeeJSCLwPPAP8jgfqrpBRCYDG4EsYIiqRvQTyPM5P91EpB1Od8Eu4EGIyvPTBRgArBOR1f51f8Hlz4/d+m+MMREiWrtcjDEm4lhBN8aYCGEF3RhjIoQVdGOMiRBW0I0xJkJYQTfGmAhhBd0YYyKEFXRjikBE4v3zfHfyOhcAEfmPiLzudR4mtFhBN6Zo/gysVNUVOStE5AMRURH59MydReRm/7asM/afm9fB/fv+vhj5vAgMFpHGxXiPiXBW0I0phIiUBwYD/8lj889ALxGpdcb6QcBut3JS1b3APOAPbsUw4ccKuok6IjJFRL7K9X0NETktIh3yecsNQAXg6zy2bQWWAvfkOl5D4FpgbAnz6+ZvsZ+57Dpj12lAcVr1JsJZQTfR6H+mKlXVwzgPH7g5n/27AqtUNSuf7aOAgf4Z9gAG4rSeS9pC/w7nEWY5SyucqVQXnLHfMqCWiLQsYRwTYaygm2h05tzTACnk/zSmC4C9BRzvE6A6zkyDscB9OEU+L91EJOXMJfcOqpqhqgdU9QBwFHgH2AE8lMe/A5xn4xoTndPnmqiXCMSLSFVVPeGf//0KYEg++1cATuR3MFVNE5EJwAM4DwQuA8wA7sxj92X8Oh92blvzOfxInF8+l6pq+hnb0nLlZ4wVdBOVclq2DXAK9d+BA8DUfPY/jNMCL8h/gFVAQ2Csqmb+2gPzP1JVdduZK/PaV0SeAm4FLlPVI3kcKyenw4XkZqKEFXQTjXIKen3/uPIBwBWqmpbP/j8CDxd0QFXd5H8YSBfyboEXi4jcjDM08QZV3ZzPbq0BH84vEmOsoJuotBfIBh7FeebldYU8sWoW8IaINFDVPQXsdz1QXlWPlSY5EWkFTAReAH4Skdr+TT7/Bdwc3YDFhTx82EQRuyhqoo5/tMpBoC1wTe6bhfLZfxOwEKclX9B+p0tbzP06AZVwuoL251py39QkQH/yHhtvopQ9gs6YIhCRK4GPgGaqejoE8rkd+CvQLsKfzWmKwVroxhSBqn4L/A1nCGMoiAPutWJucrMWujHGRAhroRtjTISwgm6MMRHCCroxxkQIK+jGGBMhrKAbY0yEsIJujDER4v8Dq/qVRL45t54AAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "pars.update(popIII_uv)\n", "pars.update(popIII_xr)\n", " \n", "sim2 = ares.simulations.Global21cm(**pars)\n", "sim2.run()\n", "\n", "ax, zax = sim.GlobalSignature(color='k')\n", "ax, zax = sim2.GlobalSignature(color='b', ax=ax)" ] }, { "cell_type": "markdown", "id": "3ddd3d3d", "metadata": {}, "source": [ "Note that the parameter file hangs onto the parameters of each population separately. To verify a few key changes, you could do: " ] }, { "cell_type": "code", "execution_count": 8, "id": "52da88af", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "4" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "len(sim2.pf.pfs)" ] }, { "cell_type": "code", "execution_count": 9, "id": "e39a05d0", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "pop_Tmin 10000.0 300\n", "pop_fstar 0.1 0.0001\n", "pop_rad_yield 1e+42 1e+42\n" ] } ], "source": [ "for key in ['pop_Tmin', 'pop_fstar', 'pop_rad_yield']:\n", " print(key, sim2.pf.pfs[0][key], sim2.pf.pfs[2][key])" ] }, { "cell_type": "markdown", "id": "3f81d3fb", "metadata": {}, "source": [ "**NOTE:** These are very simple models for PopII and PopIII stars. For more sophisticated approaches, see [More Realistic Galaxy Populations](example_pop_galaxy) and [Including Population III Stars](example_pop_popIII). \n" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.7.6" } }, "nbformat": 4, "nbformat_minor": 5 }