{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Load the training data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "from load_shamnet4_traindata import load_shamnet4_traindata\n",
    "\n",
    "traindata = load_shamnet4_traindata('shamnet4_traindata.h5')\n",
    "smf_all_batches, shmf_all_batches, scatter_all_batches, smhm_all_batches = traindata"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Compare predicted to target SMF for a particular point in the training data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZIAAAEVCAYAAADKN2OaAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nO3deVxc9b34/9eHJRCSwDAJJCELYTD7zqJG4xIDbq1+NUKi3tb7tdeArd1rwdzbX83tbYvQ3mttv60Sc7WrmoC7bVUm7juLJmZTwmRfTYbJHtbP7485Mw7DzgBnhryfjwePcOacOfM+M+S857MrrTVCCCFEX4WZHYAQQojQJolECCFEQCSRCCGECIgkEiGEEAGRRCKEECIgkkiEEEIEJMLsAHpKKVUAOAArgNZ6jbkRCSGEgBApkSiligGH1rrcSCCpSqkcs+MSQggBKhQGJCql6rXW8T7baUCx1jrbxLCEEEIQAiUSI2n4cwFZgx2LEEKI9oI+keBuE3H6Pea/LYQQwiShkEgsne1QSnW6TwghxOAIhV5bLoyeWj78t9uIiYnRZ8+e9W6PHz+epKSkAQhNCCGGjurq6qNa64TePi8UEomT9qUSC4DW2tXRE2bNmkVVVdVAxyWEEEOKUmp3X54X9FVbWusa3KUSX1bAbkI4Qggh/AR9IjGs9xs3kg2UmhWMEEKIL4VEItFa5wM2pVSWUioPqNNal5sdlxBCiNBoIwFAa11idgxCCCHaC4kSiRBCiOAliUQIIURAJJEIIYQIiCQSIYQQAZFEIoQQIiCSSIQQQgREEokQQoiASCIRQggREEkkQgghAiKJRAghREAkkQghhAiIJBIhhBABkUQihBAiIJJIhBBCBEQSiRBCiIBIIhFCCBEQSSRCCCECIolECCFEQCSRCCGECIgkEiGEEAGRRCKEECIgkkiEEEIERBKJEEKIgEgiEUIIEZDzPpF84DjGA//czsa9LrTWZocjhBAhJ8LsAMy2aZ+LR9928MibdSTFRXPNnHFcN2c86cnxhIcps8MTQoigp0LhW7hSygLkAaO11oXdHZ+RkaGrqqp6fH7XmUbs247w8uaDvFV7lMbmVsaMHEb2rHFcO2cci2yjGRZx3hfehBBDnFKqWmud0dvnBX2JRCmVBViA1IF6DUvMMHLSJ5KTPpFTDc28vv0IL285xPOf7OfJj/YwKjqCq2Ykcs3scVwxLYERUUH/tgkhxKAJ+jui1toOoJTKxJ1QBtTIqAhumJ/EDfOTONfUwrs7jvLy5kPYtx3m+U8OMCwijMunjuHqWeNYOjOR0SOjBjokIYQIakGfSMwUHRnO0pljWTpzLM0trVTuqueVLYeo2HoY+7YjhCnImGLl6lljyZ41luTRI8wOWQghBl1ItJEAKKWKAYvWOr+7Y5OSkvTBgwe92/fffz+rV6/u9PjWhgbConpestBas+XACV7dephXtxxi+6GTAEwfO4qsWYlkzxrHvAlxhEljvRAihPS1jWRIJpLeNrbvyc/n3KebiZo61fsTPWsmw+fO7dHz9zrP8OrWw1RsPUTlrnpaWjWJo6JYOnMs2bMSuSR1DNGR4T2ORwghzDBkG9sHQ+zV1xAxZgwNtTtwPfMM+swZYjIySP7rXwA49LOfoaKHE3XBBe5Ek2ojLCbG+/xJ1hj+bXEK/7Y4hfrTjbz+2REqth7mBaOxfnhkOIunjiFrZiJLZiSSOCrarEsVQoh+N+iJRCmVA6zo5jBnT0oe/cVyyzIstywDQLe20nTgIK1nTru3tebclq2c27YN3djofoJSxH/ta4z7j39Ha83JV15hWIqNqJQpxI8YxrK0iSxLm0hDcwsfOpxUbD3Mhm2Hqdh6GID5kywsnZHIVTMSmZ0Ui1JSBSaECF1StdVDuqWFxj17aKitpaG2lqhp04jNzqb56FFqF1/mPigigmHJyURNnYolN4eRl16Kbm0FrSEsjG0HT7Jh22E2bD/Cxn0utIaxsVFcNSORJdMTufSCMdK1WAhhGqnaGmAqPJyolBSiUlLg6qu9j4dbLKQ8/zwNO2pp+LyWhh07OLd1K81ffAHAuW3b2H37vzAs1YZl6lRunzqVO+dM5fTNF/LmgXO8vv0IL248yJMf7WVYeBgX2awsme6uAksZI73AhBDBL+hLJEqpNCAL8JRESgG71rqms+cMRImkrxr37KH+yado2LGDhtpamg8dAmDy448xYtEiztTU4HzmWQ5bJ1AdFs/fTw1n4yl3w3zy6BiWTE/kiukJXJwymuHDpMFeCDFwhnyvrd4IpkTir+XECRp27CBq2nTCR47g+AsvcPiXRbS4XN5jlCWejT/5Da8cVez5eAsRZ09zKD6JedMncMW0BK6YlkBqwghpWxFC9Cup2goR4bGxxKSlebfjbryR2BtuoOXoUW/7y7naWpZfs5AVkZHs2/VPTr7yJADHRlrZMWIsf40dj/2Sm7ls+liumDqaRRckEDc80qxLEkKc5ySRBAGlFBEJCUQkJDDikkva7Bv7zXziLr+MhtpaYj//nITtnzH/+Db2T7iTFzceIPn3v4CThzk+bjLDp01lcsZcpl+azvBJE026GiHE+UaqtkKQbm1FhYXR1NLK5odKOf7e+0Tu2Yn1xFEAPh+dzEv5P+eyqWO4+L0XGTM2nujp04maNo3wuDiToxdCBCtpI/Ex1BNJZ5xHnHz8Vg2bdh/jmZYE9h47w59e/QWJZ79sfwlLHMvo229lzN13A9DgcDBs0iRUpFSNCXG+k0Ti43xNJP52HzvNu7VH+eTjWg5v3EzCF/uYcvwgR6bMQH/1ZhYnRjL5GzehIiMZdsEF7lLL9OmMvOJyomw2s8MXQgwyaWwX7SSPHkHy6BHcfnEyra1L2XrwBO/VHWX3jmNU1uyj7MwZLkm/jfSWY8w6c4SEN94i/LnnCIuJIcpmo2HnTg7/sojoGTOImjGd6JkzGZacjAqXbshCiC9JIjlPhIUp5kyIY86EOPIuT6WppZVN+1y8t2Mu7zqO8dvd9TQ0t2JpPEnqnjjmv7iVy5oOMvHQYU6//z40NwOgoqOZvPZRYjIyaDp0iOZDh4iaPp2w4cNNvkIhhFmkaksA0NDcwsd7XHzocPKB4xjVe+ppbG4FYHZCNEtjzpDedJTJzv1MuiefiIQEjj32OEdKSiAsjGEpKUTPmEH0rJlYVtxK+EgZlS9EqJE2Eh+SSALX0NzCxr3H+dBxjA93OqnZU8+ZxhbAPeL+wilWFsXDvBP7iN3voGHbds5t20bzF18wvbqKsKgojj76KOc2bSJ61izvT0RCgslXJoTojLSRiH4VFRHOhSlWLkyx8h2gqaWVLQdO8NHOY1Tuqse+7TBlZ5oAGDNyBumLFpFxq5X0MRE0h0cyDKClhYbPazlZYf/yvDNmYHvuWQDObd9OxJgxRIwZM/gXKIToN1IiEX3S2qqp++IUVbvrqdzlpHKXk73OswBERYQxb2Ic6clW0pPjWWCNYMS+nZzbsgXd1MTou+4CoO4rX6Wxro6IceOInjOb4XPmEHPhRcSkLTTz0oQ4b0nVlg9JJOY4cuIc1bvrqdpdT/XuerYcOE5Ti/vvK3l0DGmT40mbbGHh5HhmjBtF48ZPOPfpp5z9dDPnNm+mcdcu4m66iaQHitBac/AnPyF62jSi584jetZMwqJlQTAhBpIkEh+SSILDuaYWPt1/nJrd9dTsqad6t4ujpxoAGB4ZztwJcSycbGHBJAsLJltIDGum9cwZIseOpfnYMXbekuOdLZnwcKKmT2PM3XcTe/XV7nVeABUWZtblCTHkSBuJCDrRkeFkTrGSOcUKuFeb3Fd/lo/3uqjZXc8ne108/u4uGlvcSWFsbJQ7qUw6xfxJccx9+VWiT9RzbvNmzm7cxNlNG70j8M9t2sSevHyGz5vH8AUL3D/z5xE+apRp1yvE+UpKJMJUDc0tbD1wgk/2uti418Une13sOnYGAKUgNWEk8ydamD8pjnkTLcwcP4qoiHAaamtx/vkvnN24kYbaWvcqlEoxZd1TDJ83j6ZDh9ANDUROnizT7QvRQ1IiESEpKiKchZPjWTg53vtY/elGNu0/zid7XGza5+LNz4/wdM0+ACLDFTPGxTJvYhzzb1zJ3G/FYYuB5i2bOfPxx0RNneo+xxNPcmzNGsLHjCFm4UKGp6cRk55O9OzZUh0mRD+TEokIelprDhw/x6a9LjbuO86mfS4+3Xeckw3u0fbRkWHMTopj7oQ45k10/0w4fYxzH37I2Y9rOFNdQ9O+fYSNGsW0D95HhYdz6s03UdHDGT5/njTiC2GQxnYfkkiGvtZWza5jp9m07zib9h3n0/0uNu8/wdkm96DJmGHhzElyTwkzb2Ics6MaGX/qKCMvzATAccONNNTWoiIjiZ43j5jMDEZefnmbRceEON9IIvEhieT81GKMbdm07zif7nPx6f7jbDlwggZjqpeRURHMSopl3oQ45seHM+MLB6M+38zZqirObdlC7HXXMeHXv0JrzbG1axk+dx7DFy4gLCrK5CsTYnBIIvEhiUR4NLe0UnvkFJ/uP87m/e7Sy7aDXyaXUVERzJ4QS9qYKOZYwpk17wKSmk5QtzQLWltRUVEMT1vIiIsXEXvtNQxLTjb5ioQYOJJIfEgiEV1pamllx5FTfLrvOJ/uP86m/e7k4pmkMjY6goyESC47u5+ZBz8jfvtGWut2MOHB/yH2uuto3LOH0++9x4hLL2XYpEkmX40Q/UcSiQ9JJKK3mlpa+fzwST7d504sm43k4hmZPyW8gdTJCcxOSWTR5jexrHkQgMhJkxhx6SWMvPxyRixeTNiwYWZehhABkUTiQxKJ6A8NzS18duik0aDvYtO+43x++CStrZoJp49y5UkHlzh3kLxnG2HNTUx86x1iR1s4u3EjavhwoqZOlTEsIqRIIvEhiUQMlDONzWw5cIKNRlfkjXtdHDh6guQTh9gZP5GpiaNY9cqDjN21HRLHErfkSkZdeQUjLr5YFv8SQW/IJhKllAXIMzYzgQqt9ZquniOJRAwm5+lG76j8j/e62LN9J9P3bOHCw9tI++Jzopsb+WLmQih+iAWTLYw4d5pwi8XssIVoZyiPbF+ltS70bCil6pRSdJdMhBgs1hHDWDIjkSUzEgFobc1k57HTfLzHxdOOI5z44CN2H29g42MfEdt4mide/k/qJ9jgkstJvul6Ji2YJVVgIqQFdYnEKI08qrXO9XmsAMjXWqd29jwpkYhgc/JcExv3HmfTlt2ol55h4tYqLnDuAeBQbCIf5XyTyVdewkUpVlLGjJDEIkwxlEskWUopm9baYWy7AJuZAQnRW6OiI1k8dQyLp46Bm9JpadVs/3QHu198hfD33uQNJ2x75lMuPLSVRS4HZxZdQcrlF3PJ1ASSR8dIYhFBrU8lEqXUFIybudb6tf4NqdvXLgMsWuvszo5JSkrSBw8e9G7ff//9rF69ehCiE6JvtNY4jp5m128fZuwzfya8pZmj0XG8PWEe26ZlMm5RJounJXBJ6hgSRslIezEwBrWxXSlVC+RqrT/p9ZMDYFR17QSWaq1rOjtOqrZEKGs5dYqTr73GoRf+QcsH7+GKS+Bb2QUcP9fMmLMuElImcfn0BC6fmkDmFCvDImQ2Y9E/Brtqa40niSilbgE82ahGa72rj+fsiUdxJ7BOk4gQoS585EgsN96I5cYbaTl5kqYDB6iZOo2tu76AZddRXxnHy+MX8MOkBZyKT+SSC8awZHoiV81IZFyczGQsBl9fSyS3aK2f9tl+AHiqJyUUpVQOsKKbw5xa63y/5xXgTlT27l5DSiRiKGo9d47jzz3P8Zde5GxVNQCHp8zksVnX89aw8QDMmRDL0hljuXr2WGaNj5W2FdErg121tUxr/UxPtwNlJB+XJ4kopbK6SiiSSMRQ17R/P8dffInjL7zA+KJfsnesjfdfq2Tbph2UM54WwphgGc7Vs8dy3ZzxpCfHEx4mSUV0bbATyQ6g2uchG+DpVaWAhVrrqb0+ccevlWWcf73xkBXI0VqXdPYcSSTifOH5/6uU4tDP/ov6J54gLDGRQxcv5YWkdF48GkZjcysJo6K4dvY4vjpvPJlTrIRJUhEdGOxEsh5Yh7srbkdytdZ39/rE7V/HAtR3sKvcd2yJP0kk4nzU2tjIqdffwPV0OaffeRdaWxmenc2mu+7j5c0HeW37Ec41tTIuNpqvzBvPTQsmMGeCVH+JLw12IknRWu/s6/6BJolEnO+aDh3i+LPPQngEY/JWorXm4B8e4ZPUDJ49onjzsy9obGklNWEENy+cwLK0iSRZZC6w891gJ5IFXTWsd7d/oEkiEaKtc59/zs6bl0FLCyMWL2bYLbm8ZpnKc5sO8dFOJ0rB4gvGkJM+kWtmjyM6MtzskIUJBjuRvAKUdXFIrtb6ml6fuJ9IIhGivabDh3GVleNav57mI0eITEpi0tq1HI4bS3nNPp6u3sd+11ksMZHckjaR2y6cxAWJo8wOWwyiwU4krwJ1fNnA7i9ba311r0/cTySRCNE53dTEydde58RLLzHhv3+NGjaMU2+/Q/i4cVQRx1Mf7eXVrYdoatEsso3mjkXJZM8aS0S4DHwc6gZ9GnmlVApgAdBaf+y3b6H/Y4NJEokQPae1pu7aa2navYcRl13G6G/cyZk5Cymr3sffPtjDftdZxsdF8/VFydx+4WQsMbIK5FBl6nokSqmFxq/azLYRD0kkQvROs9NJ/VNPUf/Ek7QcPUrUzJmMLfgx0RddzGvbj/DH93by7o5jDI8MJyd9InddlkLy6BFmhy36WVAsbGWUUh4AKrXWv+63E/eSJBIh+qa1oYETL77Iscf/yNhVqxi5+FJaTp5ERUbyWX0jj72zk+c/OUBzayvXzx3P3VekMmdCnNlhi35idolkJZCPe86tNcB6rfXxgE/cR5JIhAiMbm0FpVBKcbjkVxx/7jmsd/5f4m+7naOt4Tz27k6e+GAPJxuaWTojke8sncqCSbLqY6gzo43kKuBuYCFQjnsix53Gvlit9Yk+nbgfSCIRov+cqa7m6B8e5vS77xIWF4f1jq9jveMOTkdG85f3d/Po2w5cZ5q4fFoCP8qexnxJKCHLjGnk7biTR7tGdaVUkdZ6Va9P3E8kkQjR/85u2sTRhx/h1OuvE3fzzSQV/RKAUw3N/OX93ax5q476M01cPWssP7x6GjPGxZocsegtM6ZIKTU2PSdQxu8KKJTuv0IMTWe3bCF85EiGJSfTUFfHyddew/q1r3FaRfDYO7tY+7aDU43N5KRN5IdXT2N8nIyYDxUyRYoPSSRCDI6jpWv44sEHCU8Yw5i77yY+N5fjzfCHN+r447u7CAuDf1ucwreuvIARUaGwsvf5LSh6bQULSSRCDJ4zNTV88eBvOFNZSeTkyST+6EfEXnM1e51n+PWrn/H8JwdIHBXFfdfN4KYFE2Tm4SDW10TSp6GqSqkfK6XWKaWWGdsLlVI7lFKvGD9T+nJeIUToiUlLY/Kf/8SkNaWERUdz9hP3ULJJ1hgeunUhz3zrEsbHRfPD9RvJeeQ9th4wrR+OGCB9rdpaClR5uvgaje/3eVZNVErdK+NIhDj/6JYWdFMTYdHRnH7/feqfeJLEe39ExKTJlNfso/if23GdbeLOS6bw/expjJTqrqAyqCUSIMUniSwErL5L7wKmtY8IIcyjwsMJi3avG9904CCn3n2Xuq/ewNH/+W9umWFlw4+uYEXmJNa+s5Ps/3mT17YfNjli0R/6mkjifX5fDmzw2z/0Gl6EEL1iuWUZqS//k7ivfpVja/+XuuuuJezNDfzy5rk8861LGBUdwTf+WMX3n/qY+tONZocrAtDXROJQSj2slHoA94j2As8OY5S7EEIQmZhIUtEvmbLuKSLHJ9HsdAKQNjmeF7+zmO8uncpLmw6S/eCbVGyV0kmoCnT23xSg2qeaaynuGYHrtdav9VuUvSRtJEIEH93aClqjwsNxPfscDbW1JHz7Hra7mvlR2Ua2HTzBrZmT+MlXZ0nbiUkGpY1EKRWrlHpEKXUMeBXIBVo9+7XWG7TWT5uZRIQQwUmFhaHC3SsvNtbtwPnYYzi+egOT6zby3D2XcPcVqayr2sv1D73NJ3tdJkcreqO3VVslQAWQgbtt5DhQo5Sa39+BCSGGrsR77yX5r39BRUezd2UeR++7j3svHsv6/EW0tGpyHn6PR99y0Noqza2hoNNEopS6RSm1zDMmxPi3zChx7NRaf6y1vg/IBH6llJKJdYQQPRaTkUHKc88y5p57OPHKK5z95BMyp1j5x3cvY+nMRH7xj21840+V0hAfArptIzG692YA6UAVYNda7/I7Jg5YaebYEV/SRiJEaGk6dIjIceMAOPHqqwxPT+epz0/yXy9tI2FUFI98LZ25E2Xdk4E2YG0kRsnjUaAYdyKJV0qtVErdZfwsMBrbpVJTCNEnniTScuIEB+9bxc7/cxM3ndtJ2d2L0FpzyyPvsa5yj8lRis70uI3EmIQxy5NYtNZrtdZrAWV0+U03EstVUs0lhOiL8NhYkp98gojRo9n3zW+RuOZ/eHFlOhelWCl8+lNWv7CF5pbW7k8kBlVv+9gdV0ot01o/43nAWI/kY6XUSq31WqOaK1sppX2P6yullAV3wz5AqvGahYGeVwgRnKKnT2dK2XqO/u7/cWztWoZ9XMNj5U9TsqGOR9/eyY4jp/j97WnExUSaHaow9HocibEWiQZ+qbXeaDy2DLANRBuJUqoU9/omLmO7GlintS7p7DnSRiLE0HD6w49o+Gw71jvuAGB91V7+49lPmRQfw+N3ZpI8eoTJEQ4tgzbXltbaMyXK60qpFqVUC5A9gA3tGUCWz7YDd08xIcQQN+KiC71J5NS777Loj8U8cessnGcaWfaH92S8SZDo0xQpWus1Wmsr7skaw7XW3+znuHxfK11rXe7zUBrusSxCiPNI04EDnHrzLeK/fxdlV8YRExXOrWvexy5Tq5iur3NtAeCZGmWwKKUKcHc/XtPVcQcOHEAp5f1ZvXr14AQohBgw8bm5TPnbX0FrWu5ZyRNjDzJt7Cjy/lJFWdVes8M7rwWUSAaLUsqilMozNuu6Oz4pKQmttfdHEokQQ8PwefNIeeZpYi66iFNFP+fx2S1cesEYfly+if99R1avMEvILbWrlKoAXFrr3M6OkcZ2IYY23dLCiZdfJvb662lsaeV7T9Tw8tYjfHfpVH6QNRWlZDnfvuhrY/ugT7GplMoBVnRzmFNrnW90/c3z66FVgXtwpBDiPKXCw4n7ylfcG3v38OPyn5O89E5+u6GWppZWCq6ZLslkEAWcSIzBh1m4e1JZcI9wr8TdltFucWaj4bzc//FOZADFSqk1nu6/QgjhSzc20nryJDf/8WdYln2L4jegVWvuu3aGJJNBElAbiVLqFtwzAqfi7pZbY/ybijsBLAvk/FprOz5jSAzZxmsKIQTR06eTsn4dUdOmcuXf/pufN26i9I06HvjndkKt6j5UBVwi0Vrf3dk+Y6GrQJUbvbUARgMVXQ1GFEKcfyISEkj+0584UFBI+j/+zE9z7+FnbymiIsP5YfY0s8Mb8ga6jSTg6Tq11g6kBCKE6EZYdDQTfvMg9X97gjtuuYVtf/+c326oZWRUOHmXp5od3pAWaCJxGFOm1AHHjMdGG//agNIAzy+EED2mwsKwfv1rAPwiewqZz63l/ufPMSIqgn+5KNnk6IaugBKJMWHjcqOtxGY87AQcWutVgQYnhBB91bhpI3M/eZ3fJ9TxA92CNWYY180db3ZYQ1K/VG1prZ/uj/MIIUR/GXn55Uz8f79Dff8H/Pb9hylUrYz+9tVcmGI1O7QhZ0BHtgfaa0sIIQIxaskSJq99lLENJ3jg7Yf58cMV1B4+aXZYQ05AJRKl1AN03qCucC/PG/CaJEII0VcxmZlMXrsW9dP7GR6h+NfHPuK5b19K4qhos0MbMgKt2vLMe+UINBAhhBgoMWkLmfbCc/z64EmWP/wu3y99g8e+l010ZLjZoQ0JgSaS9UCG1npDRztlVKkQIliosDDmTIjj8dPvcOLld7g/LooH7rpS7lP9IOBp5DtLIsb+TvcJIYQZ5n49l3FNJ7ni0Z/xyAs1ZoczJITENPJCCNFfYtIWkvLIH5h05hhJRfdhr5Ka+UD1KpEopR7uZv8DSqlXlFL3BhaWEEIMnJGLFjHhoYdIPXGQwz++l7oj0pMrEL0tkZQDKKUWGLP+eimlioCjWutrgI+l668QIphZly5hxE9X84+ZV/HNv9VwuqHZ7JBCVq+rtpRSO3AnlJ1+JY8crfWvwds2Ii1YQoiglnJrDt/8wa3sOHKKX/3ueZktuI96m0jSgHSt9QVa69HAcZ+SiX/ikPVDhBBBb/HUMfxiSiMr1v6EV37+O7PDCUm9TSQOrfVxn+31uBefAqj3O1ZSuxAiJCz/xg3UXrCQiX97hC3PvWx2OCGnt4nEqpS6SikVq5SaAtwHVCml4oB4v2Nt/k8WQohgFB4ZwaLHfs/++CQa/r9V1G/eZnZIIaVXiURr/SiwHNgFlOEe0Z4NPADkK6XuVUpNUUrdhYx2F0KEkITEeCy/+S1nwodRe1c+LadOmx1SyOj1yHZjRUT/VRGfBlBKOXGXUl7VWr8WeHhCCDF4LrpoFo9956dUvfoeh3a4uHHBCLNDCgl9miLFKHGkA1bcS9+uBe/6JJ0uvSuEEMHujv97HS82WPjJc5tJj1dMSE4yO6Sg15fuv6/ibmB3AFVAhlKq0n9ciRBChKKI8DAeXL6A5C9288UNX+H4qxVmhxT0elUiMUoiuX49t1BKWYA84Nf9GJsQQphiypgR3H7bVeyqLie8oJDpM6YzbPJks8MKWr0tkdT7JxEArbUL2Nk/IQkhhPluvTSVN277AedawPG9H6AbG80OKWj1NpF0NTZExo0IIYYMpRT/fucSSi+8Fb1tK4d/85DZIQWt3iaS0cb4kTaUUguA1P4ISAghgkVibDRXf/M2XkpZxPadR2QKlU70qo1Ea/2oUmq9UiqFL8eJ2HCPeF/R79EJIYTJctMncsctK3l8z3FeqT/LJGuM2SEFnV732tJaL8fdsG7H3WijWTMAABqiSURBVGsrT5KIEGKoUkpRdMt8FPC7PzzP4ZJfScnET58WttJaf6y1flRr/Stj7AiDtQaJUqp0MF5HCCE8JsbHUHDtDHRNFc7HHuPE3/9hdkhBpdOqLaXUK704j8I9QHFAu/8qpYqRObyEECb42sXJ3HLljew4shX1s58Rc2EmkYmJZocVFLpqI1FAIT2bDl7hnm9rwCil0gby/EII0ZXwMMV/LVvAtz5fzh/eeJBDP72fiQ//AaVk6aWuEkmhp9qqJ5RShf0QT1cygArca6IIIcSgmzsxjquyM3n80LWsfONFTr3+BqOuWmJ2WKbrtI2kN0nEOH7ABiQqpXJwr33SIwcOHEAp5f1ZvXr1QIUmhDjP/Oia6bw7fyllS+5g+KWXmh1OUOjTpI2DyZh+xaW1dvW0CJmUlMSBAwcGNjAhxHkpNjqSf79hDt97qpnpmw6xfP441LBhZodlqj712hpky7XWdrODEEIIjxvnJ5GRHE/ZExXUZmVz9tNPzQ7JVINeIjGqqbobd+LUWucbDeySRIQQQUUpxf03zGbFbw5y6lwzB++/n5SyMlR4uNmhmWLQE4nWuhwo7+HhViDLp0orE7AppQqAcq21rMIohDDF3IlxfOWiC/jdkesp+OivuJ55hvjcXLPDMkVQt5EYVVreEolSKg+waa1LzItKCCHcfnztdJZsOsD+Ax8S/puHiL3uOsJHjjQ7rEEXCm0kgDeJ5GKUSIxGeCGEME3iqGi+vXQaJanX0XLsGCde+rvZIZlCDcU5YzIyMnRVVZXZYQghzgPnmlpY8us3WHD2EL9ffRthYSHz/bwdpVS11jqjt88L3SsWQoggEB0Zzg+zp/HPRgv/3HyYFldPJgMZWiSRCCFEgJalTWT62FE8//jz1F655LzrDiyJRAghAhQepii8bjrvhCXQGDGML377O7NDGlSSSIQQoh8smZ7I3GlJlE9dwum33+ZMdbXZIQ0aSSRCCNEPlFIUXDuDdRMuojHWwhcP/dbskAaNJBIhhOgn6cnxLJo9gSemXsWZykoaamvNDmlQBPWARCGECDU/yJpGztZMZlx3FTOnTjU7nEEhJRIhhOhH8ydZuGLOBB7c0czxs03o5mazQxpwkkiEEKKffT9rGifONfPBdwrZe889Zocz4CSRCCFEP5szIY5rZ4/j3eOK02++xbnPPjM7pAEliUQIIQbAd5dO5ZlJF9M8LJpj//u/ZoczoCSRCCHEAJiVFMuF86bwasrFnPj7P2jav9/skAaMJBIhhBgg9yxJ5cnkxbQCzj//2exwBowkEiGEGCDpyVZSZ9t4+NJ/ZeSd3zA7nAEjiUQMGrvdTnp6OoWFhf12rpKSjtc4q6mpCfi1XC4XhYWFrFmzhpKSEmpqajp9Pf+4+uMaB/KcvdHV++BwOCgsLEQp1el743A4iI+PJzs7m/Ly8jbPyc3N9Z63sLCQ1NRU065zoHz7qgt4yTqLF/Y0mB3KwNFaD7mf9PR0PdSVlpaaHUKfYigtLdUFBQX98vrFxcW6uLh4wF4rKyurzXZ1dbW22WzdPq8/r7E35xyov4mevA95eXmdvjcVFRU6LS1NV1RUtHkc0NXV1W0eq6+v13l5ef0QdfBobW3VN/zubX3X9x/W++5bpVtbW80OqVNAle7DPVdKJCHI4XDgMnnNg2CIoTtWq7XPz+3o2tLS0sjKygokpAEzUJ9HT9+H7OxswF0SDITFYiE9PT2gcwQbpRTfvCIVfegQJ559ljMfVZodUr+TRBKCiouLB/T8PbkhDXQMwcBut7d7L/Lz802KpmsD+Xn09H3Iyclh3bp1bR5zOBxkZPRswT2HwwGAzWbrY6TB6+rZ49g5dxFnomKof+pJs8PpdzLXVoix2+04HA4cDgcWiwWbzeb9duj5NuhwOKioqKC4uBiLxYLdbic/P5/8/HwsFgulpaVs2LABi8VCTU0Ndrsdm81GRUUF2dnZVFZWUlxcjMPhoLy8HJvNRmVlJatWrfKer7MY/LlcLoqKisjMzASgrq6u3f41a9Zgs9lwOp3ec9ntdgoLC8nIyKC0tJSamhoKCwux2WyUlpZ6n19XV0d5ebn3WnJycjq9EXV2PR2xWCzk5OSQkpJCXl4e2dnZZGVlkZaW1uZ8paWlZGZmYrFYsFqtbfaXl5cDUFFR4Y29uzj6cs6uPo+uPsPO/iZ6+z545Ofnk56e3iapORyOHicGu91OXl5e0Jb6AhEepvj6FdN4+d0Mbn61guYvviAiIcHssPqNJBLgP1/cwtYDJ0x57VlJsdx/w+weH5+VleVNGHl5eW325ebmUlZWRk5ODgCFhYWUlpaSlZVFfn4+69ato9pYI8Fzw1i6dCk7d+7EYrFQWVmJ0+n03giys7O9N/60tDRyc3OpqKjoMgZ/S5cubXODqqysbLe/2mfdhtzcXKxWqzdmz760tDTvNfhyOp3e683KyiI1NbVdsvLo7Ho6U1ZWRklJCaWlpZSUlGCxWNiwYYP3JpqdnU11dbU3ia1cudIbr91u976PTqeT0tLSbt/Xvp6zq8+jq8+ws7+J3r4PHjabzZvYepoM1q1bR1VVlfdLzFCWmzGJZdMXs6zuLVxPP82Yu+82O6R+I1VbQ0h1dbX3P7fNZvNWFXh4vhn632x8byCem86aNWva3ChsNhtVVVW9isdut7c7f2pqqvd3zzdlXytWrKCoqKjHr+H/fJvNxpo1a9od19frKSgooK6ujvr6epYvX87SpUu9sVssFu+1paWlsWHDhg7jslqt3s+iqzj6es7O9OSaO/ub6On74G/FihWUlZUB7hJyd9VaK1asIC8vj1WrVnV53FAwIiqCq67O5J9TLubEyHizw+lXUiKBXpUIgl1hYSGZmZk4nU6cTmebfZ7qJV/Lly+nvLzc+63WcxOoq6vD5XJ5kwHg3ddTNTU1XVZrVFZWtmsQt1gs3d4gu2Kz2ToskXR1PdnZ2W1e03OsJx7Pv6WlpVRVVVFTU4PD4egwdo/OGvq7iqOv5+xMTz7Djv4mfHX3PviXSvLy8khJSaG0tBSn09lhFVhH0tLS2v29DkX/umgKi9/K5URsMqvNDqYfSSIJceXl5eTk5OByuUhPT29TLeJpS+iq8Tw9PZ20tDRv9YJHZmYmNTU1Paqi8MTgLy0trV1VlK/U1NR2VUsul6vT5NOTG43T6eyw109X19NR9ZanzcH/ujIyMrBYLN1eW2e6iqOv5/Tn+Tx68xl2prv3wZ/FYiEjI8NbuuqNodg24m9cXDQ3zk/i+Q92cM/YsyRc1LOOCMFOqrZCkM1m49ixY20eq6qqalMt4vmG7bkRdKa6urrDxvKcnBycTmebJORbZdRRDP485/Q9R3V1tXc7Ly+vXXfRdevWeas5PA3wHp3d8P23PdU0vvF3dz0dKSoqapeEPY3Hvo3ZPT1fd3H09ZzQ8efRl2vuSFfvg4fv55Sbm0thYWGbv6lg7yo+mL6xOIV/qXmOw3l5tJ49a3Y4/SLoSyRKqTwgHfCUyXOBYq113+s/Qpynm6WntxO4b0IZGRnex9LS0rzfDH2/6Vosljb14RaLxdtu4XQ6vb2kbDYbZWVlbXpc+d4YOoqhIxs2bKCoqMjbkGqxWCgvL/d+a/b0PsrMzMThcJCfn++tDsnKyqKsrMzbUyk7O9s7wjovL897Q/fsr6ysbFM/X1ZWhsPh8Db+dnU9/iwWC6tWrWL9+vXex+rq6tr0GPNcm+/5/F/XarVSWlrq7T2Vk5PTZRx9PWdnn0dnr2W32zv9m+jN++AZpe7b2J+Xl+etXvT0ynM4HN6OAr497zx/G9210QwlcybE8VjaZYQ/8wEn7Buw3PBVs0MKmHIPZgxeRiIpBixADbBSa93lqKeMjAzd24bh81F5eTlOp7PNf2JPN9uuejMJIQLzfM1e4v5tOZZZ01nwtz+aHY6XUqpaa93r+raQqNrSWsdrrZXWOr27JCJ6rqM+/j1tHBVC9N2185L40JZJZE0lzUePmh1OwIK+aksMnIKCAu+EeZ72CJfL1eveWUKI3omKCMdy442Ebalgzz/s2O641eyQAhIqVVsATsAKoLXussUwKSlJHzx40Lt9//33s3r16oEKUQgheu3g8bPc9h9Pcs31F7Pq+llmhwP0vWorFEokVYDL07iulCpTSjm11uWdPSEpKYkDBw4MWoBCCNFb4+OGM2vRfNZV7eMH2dOJjgw3O6Q+C/o2Eq11jV8PrUpg6A+DFUIMeV9LS+Ibb/+F9x9aa3YoARn0EolSKgdY0c1hTq11vnF8ltba7rPPAUiLsBAi5F08fRynTx+k4R9/h3uDc2bpnhj0RGJUSXVaLeVLKWUDKpRS8Vpr3xFN5+0YEiHE0BEWpmhadDnJf3+KHZ/t5oLpyWaH1CdBXbVlVGkV+iWRFbjHlQghRMhbcNv/IQzNh3973uxQ+iwUGtvLlVIFxu+jgYruem0JIUSoGJ8+jz1xCah33qKx+dsMiwjq7/cdCvpEYpRKSsyOQwghBoJSirAbb2bHB59TsfUwX5k33uyQei30Up8QQgwx6au+zyuX5vLkR3vMDqVPJJGEGLvdTnp6OiUl5hXSampqKCkpoby8nJKSkqCd2bWmpob09HQKCws73O6My+XyTg5ZUlLivV74cpJCpVSnn4HD4SA+Pp7s7GzKy8vbPCc3N9e7Pojdbic3N9d7Ls9Mzb09vwh94WGK5WlJ7Px4K3udZ8wOp/e01kPuJz09XQ9lxcXFuri4uN/PW1pa2u0xdXV1Oisry7udl5fX7nk9OY/WWldXV+uKigpdVlY2oNdUUFDQ6XZHfK/PE6fNZmvzWF5eXrvHPCoqKnRaWpquqKho8zigq6ur2x1rsVjanaMv5xehrfbff6orZy/Qv31lq2kxAFW6D/dcKZEIwP0ttyclC7vd3mZt7eLi4jazB/f0POBer91qtZKTk0NBQUGb6cj7i/+qgt2tMthR7Glpae2mnPe8B/0d72CdXwSfxCsuZUTzOTa+8g46yKeu8ieJRAB414rojv+N1n8VvJ6eB9xrb3hmG/Zf0tVMdru93XXm57cfLOZZA8SXw+Hodp3ynhro84vgMuLii9EqjPE7NlG9u97scHol6HttDZbdX7+j3WOjrrsW6+2303r2LHvz2t9I4m6+Gcuym2mur2f/d7/Xbn/8bbcSe/31NB08yIGC9vXy1jvvZNRVS/oUb11dnXc505qaGnJycrxTwnsWPLLZbFRWVrJq1SrvDdqz8JHL5aKiooLS0lLsdru3ft5isXS4YiK4b7CedUosFgsul4t169aRlZVFcXFxj8/j4Ttl/fr169tcQ0evnZ+fT35+vnfd8A0bNnjj8FyX0+ns9nW7YrFYyMnJISUlhby8PLKzs8nKyupwev38/HzS09PbJM+Opubvq4E+vwgu4bGxRM2ZQ/r+Wp6u2UfGlK5Lz8FEEkmIcjqd3nW0s7KySE1N9a5Kl52d7f09LS2N3NxcKioqKCkpaXNT9Czp6lmFD+hypbqOjrNYLN7X6ul5fHlW/auoqOhy+vqsrCzy8/NZt24d1dXV3tcGdxWZ5zFwL/VqtVr7vLZKWVkZJSUllJaWUlJSgsViaVN68rDZbNhsNu8KjD1RWlraZk1537j99eX8IrSNWnwJUx8ppaiyjnM3zA6ZiRwlkRiS//LnTveFDR/e5f6I+Pgu90eOH9/l/r7w/1Zqs9m863H73vBsNhue1SJtNhsrV64kPz+f5cuXB8XypjabzRtHbm5ut2uheK7b8xxPycvXihUrKCoqCmhdlYKCAgoKCrw9uJYuXUp9ffvqhhUrVlBWVuZNot1VO/kuJQzukpbvMraBnl+ENstNN7EzaRrHKhup2HqYG+YnmR1Sj0gbyRBhs9moq6ujrq4Ol8uF3W73/nhuqDk5OaxatYqysjLi4+O77AabnZ1Namqq92eg5eXlYbfbKSkp6fK1PWuPe1RWVrZrQLdYLN7SVm+5XK427SOeajSbzdZhw3deXp43ETidzn5v4xno84vgMiw5mcxbriUxfiRP1+wzO5wekxLJEOF0OklPT8dqtVJTU9NpG0dOTo63Siw9Pb3DOvfy8vJ+WbO9vLzc+1r+ampqyM3N9VaLwZfJsDevnZqa2u54l8vV53YETxuPf9wZGRkd3sQtFgsZGRne9qr+NtDnF8Gn8fPPuLf+I+49MY+jpxoYMzLK7JC6JSWSEOX/jdvhcJCXl0dOTo53yVwPT5VXRUVFm+f5JhubzcaxY8cCjqun57FYLO2SncPhaNO1uCfy8vLalRTWrVvHqlXuJWv83wv/7Y4UFRW1O8Y/4TqdTu/vubm5FBYWtrmeQAdpDvT5RfA681ElM1/4M3/PTQ2JJAJSIglJnl5JnlHNlZWVbdoDysrKKCoq8lYDeW5Aqamp2O12b6+tzMxM783R09XU0/upI3a73dsdNS0tDavVSllZmTcBZGVl9eg8nmvIzc31Jrnq6mqKi4s7LcH4vrbFYmnTvlNRUUFhYSGZmZk4HA5vO0RNTY03Prvd3ibezhqwLRYLq1atatNuUVdXR2lpKfDlyHbfTgV5eXnekpWnB5nD4fD2trLZbN7nFxUVkZ+fT1ZWFna7neLiYlwuFyUlJd5r7+35pSF+aBlx6SUAjKvdBGkzTI6mZ4J+zfa+yMjI0J4GZiGECCVaa3ZcuYThCxYw8aHfDOpr93XNdqnaEkKIIKKUYsQll9B86JDZofSYVG0JIUSQGf+fq1HDhpkdRo9JiUQIIYJMKCURkEQihBAiQJJIhBBCBEQSiRBCiIBIIhFCCBEQSSRCCCECIolECCFEQCSRCCGECIgkEiGEEAEJiUSilLIopYqVUnnGT9+WvuvA6tWr++tUYhDI5xVa5PMKOX1aSSvoJ21USlmAMq11trGdB2RrrXM7e05vJm1UShHs74H4knxeoUU+r9BifF6qt88Lhbm2HgVKfbbXA3aTYhFCCOEnFKq2cgC7UsqmlErTWru01n1bR1UIIUS/C+qqLaMtpBrIBhyABVgFrNRad7pEnFLqNBDj89BB4EAnhyd1sU8EH/m8Qot8XqFlptZ6RG+fFOxVW94l9jylEKXUOtzVXZ22kfTljRBCCNE3wV615Sl1+LacO3BXdwkhhAgCg14iUUrlACu6Ocyptc7HnTTwq8ZyGeexdFW9JYQQYnAEdRsJgFKqHkjxJA2j3WSD1jre3MiEEEJAaCSSAsCltV5jbBcDx7TWJQGetxioA1KBOs/5hfmMsUN5wGitdWEH+wtwl1atAPLZmasHn1eX+8Xg6urz8NkHkAlU9OT/V7A3tqO1LlFKFRg3D+ifJFIBFGut7cZ2qVIqy7MtzKOUysLdOy+1k/3FQKXWutyzrZTK8WyLwdWDz6vL/WJw9eDzWOWbXJRSdcYgxS6TSdAnEnAnk/46l1LKBmR5RsobKoBCZKCj6XySeybuP3h/eX7fotYBxYAkEhN093n14PMUg6irz8Mojdj8nlKK+97YZSIJ9l5bA6GjebocQNZgByJ6p5M51lzIZydEf8kyvmx7uGifXNoJiRJJP6uBdr2+rB08JoKPFXD6Pea/LYToA+Pe59+JKZse1NScdyUSY2BjDW2zbL/NJiwGVKdVI0axXAjRT4z/U1m4q7a6dD6WSACWAquUUhm4v9F2NF5FBB8XRunRh/+2EKJ/PArkaq1rujvwvEwkRsLw7ZmQhVHlJYKak/alEgvIlwAh+pPRS7a0pz1Zz7uqLfD23PKVS9up6kUQMr4Z+ScMK9LbToh+Y8w+UuPTw6vbziznZSIBqj09gIx6wAwZ1BYy1ht/6B7ZyJcAIfqFkTSsQJWxMq2NHrQhB/3I9oHgcyOyAqky2jZ4GAk+C8g3HioF7L71tEax29thQr4EmKe7z6snn6cYPF19HsaX6voOnlbe1Yq0cJ4mEiGEEP3nfK3aEkII0U8kkQghhAiIJBIhhBABkUQihBAiIJJIhBBCBEQSiTgvKaXSlFIVxto0QclYa6VOKaW7mktMKWVTStUbxxbLvGNisEkiEeclYxxDUI8fMsY3eQZbdjWVt2exonytdaFMFyMGmyQScT4LhRuuiy7WhDBGIjvhy0WLhBhskkiECFLGKGQ77tmpOyuReFa1k5HiwjSSSIQIXhnG+jkOOlhj22et+h4tPiTEQDkvp5EXojNKqTy+XHXR6j+Pl89+31JArtY6n4HjwG/iPGMyPYexmYV73XohTCGJRAiDUqoMKDRKARizn1ZorbON7VKgTGttN3pGVWutUxn40kAdkOP3WJrWutxnHfuqAY5BiE5J1ZYQeNsjbJ4kAt7FshxGKQRgOcYN29MzyudG3t/x+JY42rSRGA3snuSVATikp5YwkyQSIdwy+PLG7asOSDd+d9B2aV9rJ8/pD1k+vbC8JSRj2+KTOKR9RJhOqraEcOtqEJ8neZQC+UqpItylkyL/koBxs1+O3zo3nhKPZ9toJO9RPFprh1IKwKaUsvk9NwtY2c25hBhQkkiEcLMDKzp4PBXwjH53Autxl17Wd1KdlGH865+YVnkWB1JKlSql7L2sjnLhV2oyqr8sdND110hoVt+qOiEGilRtCYF3pLvDt82jg2WYs41jO00CRnVUm5u3sSKn72N1uEsSHTJi8D+/A8j2G3SYBbg6SRY2uh4NL0S/kRKJOC8Z3+aLgQylVIHWukRrnauUKjD2gftGvNTnaWVAtVLKxZcljvIeLNVsA475bLuATKBd9ZZSqhjIA5xKKadPNZYdY7oUI758jJ5cxnOKtNYuI2nZgNGAxUhKjh5UpQnRZ7LUrhA94FnrWmtd4vd4ATDarz0kC5+xJcYxeJ5r9AJLH6ixJ0aiqcDdtpMu1VtioEmJRIieyaeDQX9a65IeziA82m97oLvr2o3XkJmAxYCTNhIheqaCDto1jNJHd4nEv0RgASr7Ka6OZGmt841SUka3RwsRIKnaEqKHjKThO1DQAu278vpXbRmPlfn22sI9gl4GEYohQRKJEP3Ik0RwlwSKPEnGeNyCUd0kjd9iKJFEIoQQIiDSRiKEECIgkkiEEEIERBKJEEKIgEgiEUIIERBJJEIIIQIiiUQIIURA/n9Jdp2l0NUXiQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "from shamnet.shamnet4_traindata import _predict_no_scatter_smf, _add_scatter_to_true_smf, _double_schechter\n",
    "from shamnet.shamnet4_traindata import LOGMH_TABLE as logmh_table, SMF_LOGSM_BINS as logsm_bins\n",
    "logsm_binmids = 0.5*(logsm_bins[:-1] + logsm_bins[1:])\n",
    "\n",
    "i = np.random.randint(0, smhm_all_batches.shape[0])\n",
    "smf_params = smf_all_batches[i, :]\n",
    "shmf_params = shmf_all_batches[i, :]\n",
    "scatter_params = scatter_all_batches[i]\n",
    "smhm_params = smhm_all_batches[i, :]\n",
    "\n",
    "logsm_bf, nd_true_bf, dlogsm_bf = _predict_no_scatter_smf(smhm_params, shmf_params, logmh_table)\n",
    "\n",
    "nd_noisy_bf = _add_scatter_to_true_smf(nd_true_bf, logsm_bf, dlogsm_bf, scatter_params, logsm_bins)\n",
    "\n",
    "fig, ax = plt.subplots(1, 1)\n",
    "xlim = ax.set_xlim(8.9, 12)\n",
    "ylim = ax.set_ylim(-6, 0)\n",
    "logsm_bin_mids = 0.5*(logsm_bins[:-1] + logsm_bins[1:])\n",
    "smf_target = _double_schechter(logsm_bin_mids, *smf_params)\n",
    "__=ax.plot(logsm_binmids, np.log10(smf_target), color=mblue, \n",
    "           label=r'${\\rm target\\ double}$-${\\rm Schechter\\ SMF}$')\n",
    "\n",
    "__=ax.plot(logsm_binmids, np.log10(nd_noisy_bf), '--', color=mred, \n",
    "           label=r'${\\rm best}$-${\\rm fit\\ 3}$-${\\rm roll\\ SMHM}$')\n",
    "\n",
    "\n",
    "xlabel = ax.set_xlabel(r'${\\log_{10} M_{\\star}}$')\n",
    "ylabel = ax.set_ylabel(r'${\\log_{10} \\phi_{\\rm SMF}}$')\n",
    "leg = ax.legend()\n",
    "\n",
    "fig.savefig('check_traindata.png', bbox_extra_artists=[xlabel, ylabel], bbox_inches='tight', dpi=200)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Grab a subset of the training data and call the batch likelihood function "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "DeviceArray(0.00532927, dtype=float32)"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from shamnet.shamnet4_traindata import log_smhm_likelihood_batch\n",
    "\n",
    "istart, iend = 100, 200\n",
    "smf_params = smf_all_batches[istart:iend, :]\n",
    "shmf_params = shmf_all_batches[istart:iend, :]\n",
    "scatter_params = scatter_all_batches[istart:iend]\n",
    "smhm_params = smhm_all_batches[istart:iend, :]\n",
    "\n",
    "log_smhm_likelihood_batch(smhm_params, smf_params, shmf_params, scatter_params, logmh_table, logsm_bins)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
