diff --git a/docs/notebooks/cone_search.ipynb b/docs/notebooks/cone_search.ipynb index 18717939..635802ba 100644 --- a/docs/notebooks/cone_search.ipynb +++ b/docs/notebooks/cone_search.ipynb @@ -7,41 +7,115 @@ "source": [ "# Cone search demo\n", "\n", - "This notebook walks through a manual, non-parallelized nearest neighbor search. This shows strategies for visualizing a catalog's partitions, and using hipscat's spatial metadata to improve performance in targeted queries." + "This notebook walks through performing a cone search of the pixels in a HiPSCat catalog. This shows strategies for visualizing a catalog's partitions, and using hipscat's spatial metadata to improve performance in targeted queries." ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 23, "metadata": {}, "outputs": [], "source": [ + "import numpy as np\n", + "from hipscat.catalog import Catalog\n", + "from hipscat import inspection\n", + "import healpy as hp\n", + "\n", "## Fill in these variables with what's relevant in your use case:\n", "\n", - "catalog_path = \"\"\n", - "ra = 24.7035278\n", - "dec = -9.3653083\n", - "radius = 2 # arcsec" + "### Change this path!!!\n", + "catalog_path = \"../../tests/data/small_sky_order1\"\n", + "\n", + "ra = 0 # degrees\n", + "dec = -80 # degrees\n", + "radius = 10 # degrees" ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 24, "metadata": {}, "outputs": [], "source": [ - "## TODO - actually implement it =]" + "## Load catalog\n", + "\n", + "catalog = Catalog.read_from_hipscat(catalog_path)" ] }, { - "attachments": {}, - "cell_type": "markdown", - "metadata": {}, + "cell_type": "code", + "execution_count": 25, + "outputs": [ + { + "data": { + "text/plain": "
", + "image/png": "iVBORw0KGgoAAAANSUhEUgAAA0QAAAIECAYAAAA5Nu72AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/P9b71AAAACXBIWXMAAA9hAAAPYQGoP6dpAAA3pElEQVR4nO3dd7gcdb348c+ek+QkpAGmAWJCgkoJAheQEjVcyCW0EFBBQEICgiihCvLjUTHEApeiSJFmAa8oUoQICKJUpd2rNAWlExSUQCAJBBMCOfP7I5w9/eSULbMzr9fz8LCZndmdnZkD33e+u3sKSZIkAQAAkEN11d4BAACAahFEAABAbgkiAAAgtwQRAACQW4IIAADILUEEAADkliACAABySxABAAC5JYgAAIDcEkRAj8yfPz8KhUJcfvnl1d6VTlViH8eNGxezZs3q1bZpOoaFQiFOPfXUau8GfXT55ZdHoVCI+fPnF5ftuOOOseOOO5b0eZqu3bPPPrukj1sNrn2giSCClHv22WfjiCOOiPHjx8fAgQNj2LBhMWnSpDj33HNj2bJlPX68Cy+8MBUDcdLpvvvui1NPPTUWL15c7V2BmvTtb3879tprrxg9erToghrRr9o7AHTu17/+dey7777R0NAQBx98cEycODFWrFgR99xzT3z5y1+Oxx9/PC699NIePeaFF14YI0aM6PXsRi0YO3ZsLFu2LPr371/tXUm9ZcuWRb9+zf8ruO+++2Lu3Lkxa9asWHPNNau3Y1Cjvva1r8WYMWNiyy23jFtvvbXauwN0gyCClHr++edj//33j7Fjx8Ydd9wR66yzTvG+2bNnxzPPPBO//vWvq7iH6VUoFGLgwIHV3o2a4DiRVW+99VYMHjy4Is+1fPnyGDBgQNTV1cXzzz8f48aNi4ULF8bIkSMr8vxA33jLHKTUmWeeGUuXLo0f/ehHrWKoyYYbbhjHHnts8c+XXXZZ7LTTTjFq1KhoaGiITTbZJC666KJW24wbNy4ef/zxuPvuu6NQKEShUCh+xuD111+PE088MTbbbLMYMmRIDBs2LHbbbbd49NFHu7W/d9xxR3z84x+PwYMHx5prrhnTp0+Pv/3tb+3Wu+uuu2LrrbeOgQMHxoQJE+KSSy6JU089NQqFwmqfY8cdd4yJEyfGgw8+GDvssEMMGjQoNthgg7j44otbrdf2MzqvvPJKjBw5MnbcccdIkqS43jPPPBODBw+Oz3zmM8Vlb7/9dsyZMyc23HDDaGhoiPXXXz9OOumkePvtt7t1HNpavHhxzJo1K4YPHx5rrrlmzJw5s9O3oz3xxBPx6U9/OtZee+0YOHBgbL311nHDDTe0WqfpsyL33ntvfOlLX4qRI0fG4MGDY5999olXX3211bp/+tOfYurUqTFixIjisTr00ENbrdPyLT2nnnpqfPnLX46IiA022KB4jcyfPz8mT54cm2++eYf7/eEPfzimTp3a5XEYN25c7LnnnsXzP2jQoNhss83irrvuioiI6667LjbbbLMYOHBgbLXVVvHwww+32v7Pf/5zzJo1q/jW0TFjxsShhx4ar732Wqv1mq6lJ554Ivbbb78YNmxYvO9974tjjz02li9f3uU+9tTLL78chxxySLz//e+PhoaGWGeddWL69OmtPsdTqdddKt25ZtpKkiQ+//nPx4ABA+K6667r87XS1jXXXBNbbbVVDBo0KEaMGBEHHXRQvPTSS63WmTVrVgwZMiSeffbZ2H333WPo0KHx2c9+NiJW/Uwff/zxMXLkyBg6dGjstdde8eKLL3b4XC+99FIceuihMXr06GhoaIhNN900fvzjH7da56677opCoRC/+MUv4mtf+1qst956scYaa8Qbb7wREavOOVBbzBBBSt14440xfvz42GGHHbq1/kUXXRSbbrpp7LXXXtGvX7+48cYb48gjj4zGxsaYPXt2RER873vfi6OPPjqGDBkSX/3qVyMiYvTo0RER8dxzz8W8efNi3333jQ022CAWLFgQl1xySUyePDn++te/xrrrrtvpc992222x2267xfjx4+PUU0+NZcuWxfnnnx+TJk2Khx56qDhAePjhh2PXXXeNddZZJ+bOnRsrV66Mb3zjGz36W9RFixbF7rvvHvvtt18ccMABcfXVV8cXv/jFGDBgQKcDt1GjRsVFF10U++67b5x//vlxzDHHRGNjY8yaNSuGDh0aF154YURENDY2xl577RX33HNPfP7zn4+NN944/vKXv8Q555wTTz31VMybN6/b+xmxaqA4ffr0uOeee+ILX/hCbLzxxnH99dfHzJkz2637+OOPx6RJk2K99daLk08+OQYPHhxXX3117L333vHLX/4y9tlnn1brH3300bHWWmvFnDlzYv78+fG9730vjjrqqLjqqqsiYlUE7rLLLjFy5Mg4+eSTY80114z58+fHdddd1+n+fvKTn4ynnnoqrrzyyjjnnHNixIgRERExcuTImDFjRhx++OHx2GOPxcSJE4vb/PGPf4ynnnoqvva1r632eDzzzDNx4IEHxhFHHBEHHXRQnH322TFt2rS4+OKL4ytf+UoceeSRERFx+umnx3777RdPPvlk1NWt+nu73/3ud/Hcc8/FIYccEmPGjCm+XfTxxx+PBx54oF1Q77fffjFu3Lg4/fTT44EHHojzzjsvFi1aFP/zP/+z2v3srk996lPx+OOPx9FHHx3jxo2LV155JX73u9/F3//+91aD4kq+7r7ozTWzcuXKOPTQQ+Oqq66K66+/PvbYY494/fXX+3ytNLn88svjkEMOiW222SZOP/30WLBgQZx77rlx7733xsMPP9zqbZ3vvvtuTJ06NT72sY/F2WefHWussUZERBx22GFxxRVXxIEHHhg77LBD3HHHHbHHHnu0e64FCxbEdtttF4VCIY466qgYOXJk3HLLLfG5z30u3njjjTjuuONarf/Nb34zBgwYECeeeGK8/fbbMWDAgG6/LiBlEiB1lixZkkREMn369G5v8+9//7vdsqlTpybjx49vtWzTTTdNJk+e3G7d5cuXJytXrmy17Pnnn08aGhqSb3zjG62WRURy2WWXFZdtscUWyahRo5LXXnutuOzRRx9N6urqkoMPPri4bNq0ackaa6yRvPTSS8VlTz/9dNKvX7+kO/85mjx5chIRyXe+853isrfffrv4/CtWrOh0H5MkSQ444IBkjTXWSJ566qnkrLPOSiIimTdvXvH+n/70p0ldXV3yhz/8odV2F198cRIRyb333ltcNnbs2GTmzJld7u+8efOSiEjOPPPM4rJ33303+fjHP95u/3beeedks802S5YvX15c1tjYmOywww7JBz/4weKyyy67LImIZMqUKUljY2Nx+fHHH5/U19cnixcvTpIkSa6//vokIpI//vGPXe5jRCRz5swp/rnpuDz//POt1lu8eHEycODA5P/9v//XavkxxxyTDB48OFm6dGmXzzN27NgkIpL77ruvuOzWW29NIiIZNGhQ8sILLxSXX3LJJUlEJHfeeWdxWUfX95VXXplERPL73/++uGzOnDlJRCR77bVXq3WPPPLIJCKSRx99tMv97K5FixYlEZGcddZZXa5XqdfddF20PG+TJ0/u8Ge9M925Zpp+ts4666zknXfeST7zmc8kgwYNSm699dbiOn29VpqsWLEiGTVqVDJx4sRk2bJlxeU33XRTEhHJ17/+9eKymTNnJhGRnHzyya0e45FHHkkiIjnyyCNbLT/wwAPbXfuf+9znknXWWSdZuHBhq3X333//ZPjw4cVzceeddyYRkYwfP77D89Pk1VdfbfccQDp5yxykUNNbL4YOHdrtbQYNGlS8vWTJkli4cGFMnjw5nnvuuViyZMlqt29oaCj+rfTKlSvjtddeiyFDhsSHP/zheOihhzrd7l//+lc88sgjMWvWrFh77bWLyz/ykY/Ef/3Xf8XNN99cfMzbbrst9t5771azTRtuuGHstttu3X6d/fr1iyOOOKL45wEDBsQRRxwRr7zySjz44INdbnvBBRfE8OHD49Of/nSccsopMWPGjJg+fXrx/muuuSY23njj2GijjWLhwoXFf3baaaeIiLjzzju7vZ8RETfffHP069cvvvjFLxaX1dfXx9FHH91qvddffz3uuOOO2G+//eLNN98sPu9rr70WU6dOjaeffrrdW4Q+//nPt5od+PjHPx4rV66MF154ISKi+DfnN910U7zzzjs92u+ODB8+PKZPnx5XXnll8W2HK1eujKuuuir23nvvbn1WY5NNNontt9+++Odtt902IiJ22mmn+MAHPtBu+XPPPVdc1vL6Xr58eSxcuDC22267iIgOr8+mWdEmTce86Xrsq0GDBsWAAQPirrvuikWLFnW5biVfd1/05JpZsWJF7LvvvnHTTTfFzTffHLvsskvxvlJcKxGr3r73yiuvxJFHHtnqs2577LFHbLTRRh1+hrLlz1pE8/k+5phjWi1vO9uTJEn88pe/jGnTpkWSJK1+/qdOnRpLlixpd7xnzpzZ6vwAtUsQQQoNGzYsIiLefPPNbm9z7733xpQpU4qf4Rk5cmR85StfiYjoVhA1NjbGOeecEx/84AejoaEhRowYESNHjow///nPXW7fNAD/8Ic/3O6+jTfeOBYuXBhvvfVWvPLKK7Fs2bLYcMMN263X0bLOrLvuuu0GVB/60IciIlp9dqMja6+9dpx33nnx5z//OYYPHx7nnXdeq/uffvrpePzxx2PkyJGt/ml6/FdeeaXb+xmx6tiss846MWTIkFbL2x6rZ555JpIkiVNOOaXdc8+ZM6fD5245kI6IWGuttSIiioPzyZMnx6c+9amYO3dujBgxIqZPnx6XXXZZrz8LFRFx8MEHx9///vf4wx/+EBGr3iq5YMGCmDFjRre2b7vPw4cPj4iI9ddfv8PlLUPj9ddfj2OPPTZGjx4dgwYNipEjR8YGG2wQER1f3x/84Adb/XnChAlRV1fX5TWydOnSePnll4v/tP1MVksNDQ1xxhlnxC233BKjR4+OT3ziE3HmmWfGyy+/XNXX3Rc9uWZOP/30mDdvXlx77bUd/q6jvl4rEV3/t2WjjTYq3t+kX79+8f73v7/dY9TV1cWECRNaLW/7mK+++mosXrw4Lr300nY/g4ccckhEtP8ZbDoPQO3zGSJIoWHDhsW6664bjz32WLfWf/bZZ2PnnXeOjTbaKL773e/G+uuvHwMGDIibb745zjnnnGhsbFztY5x22mlxyimnxKGHHhrf/OY3Y+211466uro47rjjurV9rWj6GtxFixbFiy++2OozCI2NjbHZZpvFd7/73Q63bTuALZWm43viiSd2+oHzttFYX1/f4XpNfyNfKBTi2muvjQceeCBuvPHGuPXWW+PQQw+N73znO/HAAw+0i7TumDp1aowePTquuOKK+MQnPhFXXHFFjBkzJqZMmdKt7Tvb59W9lohVnwm677774stf/nJsscUWMWTIkGhsbIxdd921W9dndz5rc/bZZ8fcuXOLfx47dmyXAXXcccfFtGnTYt68eXHrrbfGKaecEqeffnrccccdseWWW6729VXidfdET66ZqVOnxm9+85s488wzY8cdd2z3bYV9vVZ6o+Usd081HcuDDjqow8/4Raya9W7J7BBkhyCClNpzzz3j0ksvjfvvv7/V2206cuONN8bbb78dN9xwQ6u/je7oLV6dDQyvvfba+M///M/40Y9+1Gr54sWLix+u78jYsWMjIuLJJ59sd98TTzwRI0aMiMGDB8fAgQNj4MCB8cwzz7Rbr6NlnfnnP//Z7ut0n3rqqYhY/bc7/eY3v4kf/vCHcdJJJ8XPfvazmDlzZvzv//5v8ffwTJgwIR599NHYeeedS/Jh9bFjx8btt98eS5cubTWYbHusxo8fHxER/fv3L/mAcbvttovtttsuvv3tb8fPf/7z+OxnPxu/+MUv4rDDDutw/a5ed319fRx44IFx+eWXxxlnnBHz5s2Lww8/vNOBfaksWrQobr/99pg7d258/etfLy5/+umnO93m6aefbvU3+M8880w0NjZ2eY0cfPDB8bGPfaz45+4MeCdMmBAnnHBCnHDCCfH000/HFltsEd/5znfiiiuuWO22q9Ob110K3blmtttuu/jCF74Qe+65Z+y7775x/fXXt/p9VqW4Vlr+t6XpbatNnnzyyeL9q3uMxsbGePbZZ1vNCrX9GWz6BrqVK1eWNdqAdPKWOUipk046KQYPHhyHHXZYLFiwoN39zz77bJx77rkR0fw3zS3/ZnnJkiVx2WWXtdtu8ODBHX7tc319favtI1Z9pqbtZ1faWmeddWKLLbaIn/zkJ60e97HHHovf/va3sfvuuxcff8qUKTFv3rz45z//WVzvmWeeiVtuuaXL52jp3XffjUsuuaT45xUrVsQll1wSI0eOjK222qrT7RYvXhyHHXZYfPSjH43TTjstfvjDH8ZDDz0Up512WnGd/fbbL1566aX4wQ9+0G77ZcuWxVtvvdXt/YyI2H333ePdd99t9fXnK1eujPPPP7/VeqNGjYodd9wxLrnkkvjXv/7V7nG6eutWZxYtWtTufG6xxRYREV2+ba4pNDv7avAZM2bEokWL4ogjjoilS5fGQQcd1ON966mOru+IVd+a2Jnvf//7rf7cdMy7+rza+PHjY8qUKcV/Jk2a1Om6//73v9t9jfeECRNi6NChfXpbYku9ed190dNrZsqUKfGLX/wifvOb38SMGTPazVj19VrZeuutY9SoUXHxxRe3ev5bbrkl/va3v3X4TXFtNZ3vtm+PbXsM6+vr41Of+lT88pe/7HBmvjc/g0DtMEMEKTVhwoT4+c9/Hp/5zGdi4403joMPPjgmTpwYK1asiPvuuy+uueaamDVrVkRE7LLLLjFgwICYNm1acfDxgx/8IEaNGtVugL3VVlvFRRddFN/61rdiww03jFGjRsVOO+0Ue+65Z3zjG9+IQw45JHbYYYf4y1/+Ej/72c+KsxddOeuss2K33XaL7bffPj73uc8Vv3Z7+PDhxd9xE7Hqd8T89re/jUmTJsUXv/jFWLlyZVxwwQUxceLEeOSRR7p1XNZdd90444wzYv78+fGhD30orrrqqnjkkUfi0ksvjf79+3e63bHHHhuvvfZa3HbbbVFfXx+77rprHHbYYfGtb30rpk+fHptvvnnMmDEjrr766vjCF74Qd955Z0yaNClWrlwZTzzxRFx99dVx6623xtZbb92t/YyImDZtWkyaNClOPvnkmD9/fmyyySZx3XXXdfjZj+9///vxsY99LDbbbLM4/PDDY/z48bFgwYK4//7748UXX+z274Nq8pOf/CQuvPDC2GeffWLChAnx5ptvxg9+8IMYNmxYMVI70hSVX/3qV2P//feP/v37x7Rp04qhtOWWW8bEiROLX0DxH//xHz3ar94YNmxY8TM677zzTqy33nrx29/+Np5//vlOt3n++edjr732il133TXuv//+4tcud/b7cXrqqaeeip133jn222+/2GSTTaJfv35x/fXXx4IFC2L//fcvyXP05nX3RW+umb333jsuu+yyOPjgg2PYsGGt/rKir9dK//7944wzzohDDjkkJk+eHAcccEDxa7fHjRsXxx9//GofY4sttogDDjggLrzwwliyZEnssMMOcfvtt3c4K/3f//3fceedd8a2224bhx9+eGyyySbx+uuvx0MPPRS33XZbvP76693a75/+9KfxwgsvxL///e+IiPj9738f3/rWtyJiVSR2Z2YLqLBqfLUd0H1PPfVUcvjhhyfjxo1LBgwYkAwdOjSZNGlScv7557f6iuYbbrgh+chHPpIMHDgwGTduXHLGGWckP/7xj9t9Fe/LL7+c7LHHHsnQoUOTiCh+Le/y5cuTE044IVlnnXWSQYMGJZMmTUruv//+dl/d29lXWt92223JpEmTkkGDBiXDhg1Lpk2blvz1r39t93puv/32ZMstt0wGDBiQTJgwIfnhD3+YnHDCCcnAgQNXeywmT56cbLrppsmf/vSnZPvtt08GDhyYjB07Nrngggtardd2H3/1q1+1+7ruJEmSN954Ixk7dmyy+eabF7+ye8WKFckZZ5yRbLrppklDQ0Oy1lprJVtttVUyd+7cZMmSJcVtu/O120mSJK+99loyY8aMZNiwYcnw4cOTGTNmJA8//HCHx/DZZ59NDj744GTMmDFJ//79k/XWWy/Zc889k2uvvba4TtPXK7f9auSmrwJu+srmhx56KDnggAOSD3zgA0lDQ0MyatSoZM8990z+9Kc/tdouOvha4G9+85vJeuutl9TV1XX4FdxnnnlmEhHJaaedttrX32Ts2LHJHnvs0W55RCSzZ89utazlVzs3efHFF5N99tknWXPNNZPhw4cn++67b/LPf/6z3f43fe32X//61+TTn/50MnTo0GSttdZKjjrqqFZf3dxXCxcuTGbPnp1stNFGyeDBg5Phw4cn2267bXL11VdX5XWX4mu3u3PNdLSPSZIkF154YRIRyYknnthqeW+ulbauuuqqZMstt0waGhqStddeO/nsZz+bvPjii63WmTlzZjJ48OAOt1+2bFlyzDHHJO973/uSwYMHJ9OmTUv+8Y9/dHjtL1iwIJk9e3ay/vrrJ/3790/GjBmT7Lzzzsmll15aXKfpZ+2aa67p8Pmafj1AR/+0/Ep1ID0KSdJmfhygwvbee+94/PHHV/vZiB133DEWLlzY7S+boDzOPffcOP7442P+/PntvkGt2k499dSYO3duvPrqq11+9o3KSPO1AtDEZ4iAilq2bFmrPz/99NNx8803d/jVvaRPkiTxox/9KCZPnmyAS5dcK0Ct8BkioKLGjx8fs2bNivHjx8cLL7wQF110UQwYMCBOOumkau8aXXjrrbfihhtuiDvvvDP+8pe/xK9+9atq7xI98Oqrr8bKlSs7vX/AgAGtfrFyX3TnWnn99ddjxYoVnT5GfX19jBw5siT7A7A6ggioqF133TWuvPLKePnll6OhoSG23377OO2009r9Ik3S5dVXX40DDzww1lxzzfjKV74Se+21V7V3iR7YZptt2v0i05YmT54cd911V0meqzvXyic/+cm4++67O32M1f0OKIBS8hkiAMi4e++9t93bVVtaa621uvza+lJ78MEHY9GiRZ3eP2jQoC6/9hyglAQRAACQW75UAQAAyC1BBAAA5JYgAgAAcksQAQAAuSWIAACA3PJ7iAAybO7cudXehUyYM2dOtXcBgDLxtdsAKSZoskFQAaSXIAKoEHFDT4gogMoQRAB9IHJIA/EE0HuCCKATYocsEU0AHRNEQG4JHmgmmIC8EkRAJokdKD3RBGSRIAJqluiB9BBLQK0SRECqiR6ofWIJSDNBBKSC8IH8EUpAGggioKKED7A6QgmoJEEElI34AUpFJAHlIoiAPhM+QLUIJaCvBBHQI+IHSDuRBPSEIAI6JX6ArBBJQGcEERAR4gfIH5EERAgiyC0BBNCaQIJ8EkSQA+IHoHdEEmSfIIKMET8A5SWSIFsEEdQ4AQRQXQIJapsgghojgADSTSBBbRFEkHICCKC2CSRIN0EEKSOAALJNIEG6CCKoMgEEkG8CCapLEEGFCSAAuiKQoLIEEVSACAKgN8QRlJ8ggjIQQACUg0CC0hNEUCIiCIBKEkdQGoII+kAEAZAG4gh6TxBBDwggAGqBQILuE0SwGiIIgFomjqBrggg6IIIAyCJxBO0JIniPCAIgT8QRrFJX7R0AAACoFjNE5JpZIQAwW0S+CSJyRwQBQOfEEXkjiMgFEQQAPSeOyANBRGaJIAAoHXFEVgkiMkMAAUDlCCSyQhBR84QQAFSPMKLWCSJqkggCgPQRR9QiQUTNEEEAUDvEEbVCEJF6QggAapcwIu0EEakkggAge8QRaSSISBUhBADZJ4xIE0FE1YkgAMgvcUS1CSKqRggBAE2EEdUiiKgoEQQArI44opIEERUhhACAnhJGVIIgoqyEEADQV8KIchJElJwIAgDKRRxRaoKIkhFCAEClCCNKRRDRZ0IIAKgWYURfCSJ6TQgBAGkhjOgtQUSPCSEAIK2EET0liOg2IQQA1AphRHcJIrokggCAWieO6IogokNCCADIGmFERwQRrQghACDrhBEtCSIiQggBAPkjjIgQRLknhACAvBNG+SaIckoIAQC0Jozyqa7aO0DliSEAgPaMkfLJDFGO+CEHAOges0X5IYhyQAgBAPSOMMo+QZRhQggAoDSEUXYJogwSQgAA5SGMskcQZYgQAgCoDGGUHYIoA4QQAEDliaJsEEQ1TAgBAKSDOKpdgqgGCSEAgHQSRrVHENUQIQQAUBuEUe2oq/YO0D1iCACgdhi71Q4zRCnnhwkAoLaZLUo3QZRSQggAIFuEUTp5y1wKiSEAgOwxxksnM0Qp4ocEACAfzBalhxmilBBDAAD5YeyXHmaIqswPAwBAvpktqi4zRFUkhgAAMCasLjNEVeCiBwCgI2aLKs8MUYWJIQAAOmOsWHlmiCrExQ0AQE+YLaoMM0QVIIYAAOgpY8jKMENURi5iAABKwWxR+ZghKhMxBABAqRhblo8ZohJzsQIAUE5mi0rLDFEJiSEAAMrNmLO0BFGJuDABAKgUY8/S8Za5PnIxAgBQTd5C1zdmiPpADAEAUG3GpH0jiHrJhQcAQFoYm/aet8z1kIsNAIA08xa6njFD1ANiCACAtDNm7RlB1E0uLAAAaoWxa/cJom5wQQEAUGuMYbvHZ4i64CICACALfK6oc2aIOiGGAADICmPbzgmiDrhgAADIGmPcjgmiNlwoAABklbFue4IIAADILUHUgmIGACDrjHlbE0TvcWEAAJAXxr7Ncv+12y4GAADyLO9fyZ3rGSIxBABA3uV9TJzbIMr7iQcAgCZ5HhvnMojyfMIBAKAjeR0j5y6I8nqiAQBgdfI4Vs5VEOXxBAMAQE/kbcycmyDK24kFAIDeytPYORdBlKcTCgAApZCXMXTmgygvJxIAAEotD2PpTAdRHk4gAACUU9bH1JkNoqyfOAAAqJQsj60zGURZPmEAAFANWR1jZy6IsnqiAACg2rI41s5UEGXxBAEAQJpkbcydmSDK2okBAIC0ytLYOxNBlKUTAgAAtSArY/CaD6KsnAgAAKg1WRiL13QQZeEEAABALav1MXlNBxEAAEBfFJIkSaq9Ez1V6xUKAABZNGfOnGrvQo/V3AyRGAIAgHSqxbF6TQVRLR5gAADIk1obs9dMENXagQUAgLyqpbF7TQRRLR1QAACgdsbwNRFEAAAA5ZD6IKqVsgQAAFqrhbF8qoOoFg4gAADQubSP6VMbRGk/cAAAQPekeWyf2iACAAAot1QGUZoLEgAA6Lm0jvFTF0RpPVAAAEDfpHGsn6ogSuMBAgAASidtY/5UBREAAEAlpSaI0laKAABAeaRp7J+KIErTAQEAAMovLQ2QiiACAACohqoHUVrKEAAAqKw0tEBVgygNBwAAAKieajdB1WeIAAAAqqVqQVTtEgQAANKhmm1QlSASQwAAQEvVagRvmQMAAHKr4kFkdggAAOhINVrBDBEAAJBbFQ0is0MAAEBXKt0MZogAAIDcqlgQmR0CAAC6o5LtYIYIAADIrYoEkdkhAACgJyrVEGUPIjEEAAD0RiVawlvmAACA3BJEAABAbpU1iLxdDgAA6ItyN4UZIgAAILfKFkRmhwAAgFIoZ1uYIQIAAHKrLEFkdggAACilcjWGGSIAACC3Sh5EZocAAIByKEdrmCECAABySxABAAC5VdIg8nY5AACgnErdHGaIAACA3BJEAABAbpUsiLxdDgAAqIRStocZIgAAILcEEQAAkFuCCAAAyK2SBJHPDwEAAJVUqgYxQwQAAOSWIAIAAHKrz0Hk7XIAAEA1lKJFzBABAAC5JYgAAIDcEkQAAEBu9SmIfH4IAACopr42iRkiAAAgtwQRAACQW4IIAADIrV4Hkc8PAQAAadCXNjFDBAAA5JYgAgAAcksQAQAAuSWIAACA3BJEAABAbvUqiHzDHAAAkCa9bRQzRAAAQG4JIgAAILcEEQAAkFuCCAAAyC1BBAAA5JYgAgAAcksQAQAAuSWIAACA3BJEAABAbgkiAAAgtwQRAACQW4IIAADILUEEAADkliACAIouX3NItXcBoKIEEQAQEc0xdPmaQ4QRkBuCCADoMIBEEZAHgggA6JQoArJOEAFAzq0uekQRkGWCCAByrLuxI4qArBJEAJBTPY0cUQRkkSACgBzqbdyIIiBrBBEA5Exfo0YUAVkiiAAgR0oVM6IIyApBBAA5UeqIEUVAFggiAMiBcsWLKAJqnSACgIwrd7SIIqCWCSIAyLBKxYooAmpVr4Jozpw5pd4PAKDERAqQJ71tFDNEAJBB1YghAQbUIkEEABlTzTARRUCtEUQAkCFpCJI07ANAdwkiAMiINIVImvYFoCuCCAAyII0BksZ9Amir10Hkm+YAIB3SHB5p3jcgO/rSJmaIAKCG1UJw1MI+AvkliACgRtVSaNTSvgL5IogAoAYJDIDS6FMQ+RwRAFRercZQre43kG59bRIzRABQQ2o9Kmp9/4HsEUQAUCOyEhNZeR1ANggiAKgBWYuIrL0eoHb1OYh8jggAyiur8ZDV1wVUTilaxAwRAKRY1qMh668PSD9BBAApJRYAyq8kQeRtcwBQWnmKoTy9VqB0StUgZogAIGXyGAh5fM1AOggiAEiRPIdBnl87UD2CCAAAyK2SBZHPEQFA35ghcQyA7ille5ghAoAUEALNHAugkgQRAFSZAAConpIGkbfNAUDPiKGOOS5AZ0rdHGaIAKBKDPq75vgAlSCIAKAKDPa7x3ECyq3kQeRtcwDQNYN8gN4pR2uYIQKAChJDPeeYAeVUliAySwQA7RnY955jB5SrMcwQAUAFGND3nWMIlEPZgsgsEQCsYiAP0DflbAszRABQRmKotBxPoNTKGkRmiQDIM4P38nBcIV/K3RRmiACgDAzay8vxBUpFEAFAiRmsA9SOsgeRt80BkCdiqHIca8i+SrRERWaIRBEAeWCAXnmOOWRXpRrCW+YAoAQMzKvHsQf6omJBZJYIgKwyIAcorUq2gxkiAKDmiVKgtyoaRGaJAMgaA/H0cC4gGyrdDGaIAKCXDMABal/Fg8gsEQBZIIbSyXmB2laNVjBDBAA9ZNCdbs4P0BNVCSKzRADUKoNtgPKoViNUbYZIFAFQa8RQ7XCuoLZUsw28ZQ4AusEAGyCbqhpEZokAqAViqDY5b1Abqt0EVZ8hqvYBAICuGFTXNucP0i0NLVD1IAKAtDKYBsi+VARRGsoQAFoSQ9nhXEI6paUBUhFEEek5IABgAA1QXmka+6cmiAAgDcRQNjmvQGdSFURpKkUA8segOducX0iHtI35UxVEEek7QAAAQGmkcaxfSJIkqfZOdGTu3LnV3gUAcsLMQb7MWry02rsAuZTGGIpI4QwRAFSSGALIt9QGUVoLEoDsEEP55LxD5aV5bJ/aIIpI94EDoLYZFANURtrH9KkOooj0H0AAao8YwjUAlVELY/nUBxEAlJKBMAAt1UQQ1UJZApB+YoiWXA9QXrUyhq+JIIqonQMKQDoZ/AJUTi2N3WsmiCJq68ACAOknlKH0am3MXlNBFFF7BxiA6jPoBaiMWhyrF5IkSaq9E701d+7cau8CACknhlidWYuXVnsXoObVYgg1qbkZIgDoLjFEd7hOIN9qOohquUQBKC+DXIDKqPUxeU0HUUTtnwAASk8M0VOuGeidLIzFaz6IIrJxIgAoDQNbgMrIyhg8E0EUkZ0TAgAAaZelsXdmgigiWycGgJ4zO0RfuH6ge7I25s5UEEVk7wQB0D0GswDll8WxduaCKCKbJwqAzokhgPLL6hg7k0EUkd0TBkBrYohScj1Bx7I8ts5sEEVk+8QBYPAKUAlZH1NnOogisn8CAQCgXPIwls58EEXk40QC5I3ZIcrFtQWr5GUMnYsgisjPCQXIAwNWgPLK09g5N0EUka8TC5BVYgigvPI2Zs5VEEXk7wQDZIkYAiivPI6VcxdEEfk80QBA94lv8iivY+RcBlFEfk84QK0yQAUonzyPjXMbRBH5PvEAtUQMAZRP3sfEhSRJkmrvRBrMnTu32rsAQAfEENUya/HSau8ClFXeQ6hJrmeIWnJBAACQF8a+zQQRAACQW4KoBaUMkC7eLkc1uf7IKmPe1gRRGy4QgHQwGAUoPWPd9gRRB1woANUlhgBKzxi3Y4KoEy4YAACywti2c752uxt8JTdA5ZgdIk189Ta1TgitnhmibnAhAVSGGAIoHWPY7hFE3eSCAgCgVhi7dp8g6gEXFkD5mB0CKA1j1p7xGaJe8rkigNIRQ6SVzxBRS4RQ75gh6iUXHAAAaWFs2nuCqA9ceAB9Z3YIoG+MSfvGW+ZKxFvoAHpHEJFm3jJHmgmh0jBDVCIuSADIFjFEmhl7lo4gKiEXJgBkgxgizYw5S8tb5srEW+gAusdb5kgjQUQaCaHyMENUJi5YAKhNYog0MrYsHzNEFWC2CKBzZohIEzFE2gih8jNDVAEuZABIPzFE2hhDVoYZogozWwTQmhki0kIQkRZCqLLMEFWYCxwA0kcMkRbGipVnhqiKzBYBmCGi+sQQaSCEqscMURW58AGgusQQaWBMWF1miFLCbBGQV2aIqCZBRDUJoXQwQ5QSfiAAoLLEENVk7JceZohSyGwRkCdmiKgGMUS1CKH0MUOUQn5QAKB8xBDVYoyXTmaIUs5sEZB1ZoioJDFENQihdBNENUIYAVkliKgkQUQlCaHa4C1zNcIPFAD0jRiikozdaocZohpktgjIEjNEVIIYolKEUO0RRDVMGAFZIIgoNzFEJQih2iWIMkAYAbVMEFFOYohyEkHZIIgyRBgBtUgQUS5iiHISQ9khiDJIGAG1RBBRDmKIchFC2SOIMkwYAbVAEFEOgohSE0LZJYhyQBgBaSaIKDUxRCkJoewTRDkijIA0EkSUkhiiVIRQfgiiHBJGQJoIIkpFDFEKQih/6qq9A1SeH3QAskYMUQrGSPlkhijnzBYB1WaGiL4SQ/SVEMo3QURECCOgegQRfSGG6AshRIQgog1hBFSaIKK3xBC9JYRoSRDRIWEEVIogojfEEL0hhOiIIKJLwggoN0FET4khekoI0RVBRLeJI6AcBBE9IYboLhFEdwkiekwYAaUkiOguMUR3CCF6ShDRa8IIKAVBRHeIIVZHCNFbgog+E0ZAXwgiuiKEWB0hRF8JIkpGGAG9IYjojBiiK0KIUhFElJwwAnpCENERMURnhBClJogoK3EErI4goi0xRFsiiHISRFSEMAI6I4hoSQzRkhCiEgQRFSWMgLYEEU3EEE2EEJUkiKgacQRECCJWEUOIIKpFEFF1wgjyTRDlmxBCCFFtgohUEUeQP4Iov8RQfokg0kQQkUrCCPJDEOWTGMonIUQaCSJSTxxBtgmi/BFD+SKCSDtBRM0QRpBNgig/hFC+CCFqhSCiJokjyA5BlA9iKB9EELVIEFHzxBHUNkGUfWIo20QQtU4QkRnCCGqTIMouIZRtQoisEERklkCC2iCIskkMZY8AIqsEEbkgjiC9BFG2CKFsEUHkgSAid8QRpIsgyg4xlA0iiLwRROSaOILqE0TZIIZqmwgizwQRvEccQXUIotomhGqXCIJV6qq9AwAAANVihgg6YLYIKscMUe0yO1R7zApBe4IIVkMcQXkJotojhGqLCIKuCSLoAXEEpSeIaocQqh0iCLpPEEEfCCToO0FUG8RQugkg6D1BBCUijqB3BFG6CaH0EkFQGoIIykAcQfcJovQRQeklgqD0BBFUgECCzgmidBFD6SKAoPwEEVSYOILWBFE6CKH0EEFQWYIIqkwgkXeCqLqEUPUJIKguQQQpI5DIG0FUHUKoegQQpIsggpQTSGSdIKosIVR5AgjSTRBBjRFIZI0gqgwhVDkCCGqLIIIaJ5CodYKovIRQ+QkgqG2CCDJGIFFrBFF5CKHyEUCQLYIIckAkkWaCqLSEUGmJH8g+QQQ5JZJIC0FUGkKo78QP5JMgAiJCIFE9gqj3RFDfCCAgQhABXRBJVIIg6jkh1HPiB+iMIAJ6RCRRaoKo+4RQ94gfoCcEEdBnIom+EERdE0FdEz9AXwkioGyEEt0hiNoTQe0JH6BcBBFQUSKJtgTRKiKomfgBKkkQAakglPIrz0GU9wgSPkAaCCIg1YRS9uUtiPIYQcIHSDNBBNQssZQNeQiiPESQ6AFqlSACMkks1Y4sBlFWA0j0AFkkiIDcEk3pkJUgqvUIEjtAXgkigE4IpsqoxSCqxfgRPAAdE0QAfSCa+q4WgijtASR2AHpPEAFUiHjqWNqCKC3xI3IAKkMQAaRYHiKq0kFUjeARNwDpJYgAMqwWgqrUQVSO4BE0ANkliAAAgNyqq/YOAAAAVIsgAgAAcksQAQAAuSWIAACA3BJEAABAbgkiAAAgtwQRAACQW4IIAADILUEEAADkliACAABySxABAAC5JYgAAIDcEkQAAEBuCSIAACC3BBEAAJBbgggAAMgtQQQAAOSWIAIAAHJLEAEAALkliAAAgNwSRAAAQG4JIgAAILcEEQAAkFuCCAAAyC1BBAAA5Fa/au8AZNHy5ctjxYoV1d4NADJmwIABMXDgwGrvBmSKIIISW758eQwftFasiOXV3hUAMmbMmDHx/PPPiyIoIUEEJbZixYpYEcvjY7F79Cs0RKGusOqOQl2L2+/9u65QvF2oq2uxvK55vab7C3Wr1m+5faHQZt1osbztum2373xfkkKh+Q21LZ+ry/vfW/7e7Vbrtl1W1/L+Fo/z3u2Wj5+0Wt68bhJNt6P4uprvb7Nui+XF9Vo+f9Pyuo63L2q1fSe3OzgWrZ6/w3U7uB2d3N9mXzrdvrNlq3n+Jp0u6+CxOjoWUUi6vS9RWPWqiq+t3bpJx8/Zcnnx+ZuXFbraPpIWPwLNz1zoZPum5a0es8X2hTbbt/hxibpW2zffX9diWd17r77l4zStW9fmdkREXbReVtfB7abH6uz+pudctayx+bmi7f2NUd9im+Z1Vz1+fSRRaLF987otlrW8/d66Tc9TX2gsPmb9e8/X9LjF5+rgseoLjcV9rG+xXtOPcX20fNymbZLiY63avmm75sepb/H661vsS9N5KT5WNB/LpvtXLYvmY9W0L4WI+vfOSPOyQtQVlzXfri80LatrsWzV7TfebIyxW82PFStWCCIoIUEEZdIv+ke/Qv8otAiSlrdX/bt51FQotAiilvfXtbh/tUFUaLddp0FUaHN/n4Oo0BwBrYKnOQxKGkQtB8ZlDqKOt+/kdnH7QovbLbZfTZCsLiLKHkQd3N+kbRB1dVx6FUQtnqviQdTydgfbtw2itts0Xy7NA+deBVEH2/Q2iFoHT/eDqOXyVf/uLIhaRkDPg6iuwyBKOrnd/SCqL+5XIereO6BNMbIqiJpuF1pEStJiWdLisaK4L83P37ysqyCq70EQ1XcriHzsG8rFTxcAAJBbgggAAMgtQQQAAOSWIAIAAHJLEAEAALkliAAAgNwSRAAAQG4JIgAAILcEEQAAkFuCCAAAyC1BBAAA5JYgAgAAcqtftXcAsurdeCciqYtCUnhvScvb7/07KRRvF5K6FsvrmtdrfO/+Qov7Cy3/Xdd8u/jwhQ7Wbbt9m/uTQvF2UihEJNFm+9Xd/97yiIjGaL1u00tuWlbX8v4Wj1PclULxZSWtljev2/T0UYiIuraP1WbdDg97i206OOyt9is62r6T222PRaH18o7X7eB2dHJ/m33pdPvOlq3m+Zt0uqyDx+roWEQh6fa+RGHVqyq+tnbrJh0/Z8vlxedvXlboavtIWvwIND9zoZPtm5a3eswW2xfabN/ixyWSVts335+0WJa89+qbHqexxTZ1bW5HRNRF62V1Hdwu/penk/vrouWyxubnirb3N0Z9i22a1131+PWRRKHF9s3rtljW8vZ76zY9T32hsfiY9e89X9PjFp+rg8eqLzQW97G+xXpNP8b10fJxm7ZJio+1avum7Zofp77F669vsS9N56r4WNF8LJvuX7Usmo9V074UIurfOyPNywpRV1zWfLu+0Lxe87JVj/fGm40BlJ4gghJLkiSGDBkS9yy9edUIb2W19wiArBgyZEgkSbL6FYFuE0RQYoVCIZYuXRr/+Mc/YtiwYdXeHQAy4o033oj1118/Ci1nroE+E0RQJsOGDRNEAAAp50sVAACA3BJEAABAbgkiKLGGhoaYM2dONDQ0VHtXAMgQ/3+B8igkvqoEAADIKTNEAABAbgkiAAAgtwQRAACQW4IIAADILUEEAADkliCCEvr+978f48aNi4EDB8a2224b//d//1ftXQKgxp1++umxzTbbxNChQ2PUqFGx9957x5NPPlnt3YLMEERQIldddVV86Utfijlz5sRDDz0Um2++eUydOjVeeeWVau8aADXs7rvvjtmzZ8cDDzwQv/vd7+Kdd96JXXbZJd56661q7xpkgt9DBCWy7bbbxjbbbBMXXHBBREQ0NjbG+uuvH0cffXScfPLJVd47ALLi1VdfjVGjRsXdd98dn/jEJ6q9O1DzzBBBCaxYsSIefPDBmDJlSnFZXV1dTJkyJe6///4q7hkAWbNkyZKIiFh77bWrvCeQDYIISmDhwoWxcuXKGD16dKvlo0ePjpdffrlKewVA1jQ2NsZxxx0XkyZNiokTJ1Z7dyAT+lV7BwAA6J7Zs2fHY489Fvfcc0+1dwUyQxBBCYwYMSLq6+tjwYIFrZYvWLAgxowZU6W9AiBLjjrqqLjpppvi97//fbz//e+v9u5AZnjLHJTAgAEDYquttorbb7+9uKyxsTFuv/322H777au4ZwDUuiRJ4qijjorrr78+7rjjjthggw2qvUuQKWaIoES+9KUvxcyZM2PrrbeOj370o/G9730v3nrrrTjkkEOqvWsA1LDZs2fHz3/+8/jVr34VQ4cOLX42dfjw4TFo0KAq7x3UPl+7DSV0wQUXxFlnnRUvv/xybLHFFnHeeefFtttuW+3dAqCGFQqFDpdfdtllMWvWrMruDGSQIAIAAHLLZ4gAAIDcEkQAAEBuCSIAACC3BBEAAJBbgggAAMgtQQQAAOSWIAIAAHJLEAEAALkliAAAgNwSRAAAQG4JIgAAILf+PyyOUkEc6PhtAAAAAElFTkSuQmCC" + }, + "metadata": {}, + "output_type": "display_data" + } + ], "source": [ - "## Refine the cone search parameters\n", + "## Plot catalog pixels\n", "\n", - "You can visualize the partitions that are found with a cone, and refine the parameters, if you'd like.\n" - ] + "inspection.plot_pixels(catalog)" + ], + "metadata": { + "collapsed": false + } + }, + { + "cell_type": "code", + "execution_count": 29, + "outputs": [ + { + "data": { + "text/plain": "
", + "image/png": "iVBORw0KGgoAAAANSUhEUgAAA0QAAAIECAYAAAA5Nu72AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/P9b71AAAACXBIWXMAAA9hAAAPYQGoP6dpAAAlJklEQVR4nO3df5DcdX348dfuwV0oIU0ov5om5CAZMYCGGhqk4FBIKCPBEUfxx4wCWiqdUalFZVBnyljbzmA7KC0ijVMRnWlHghVjYUAHgxFlZEbFQVorhKAdvyU02gQSDCe5z/ePy+3t3u1e7sfufn68H48ZJpvPfj6ffd9n9+T99L271LIsywIAACBB9bwHAAAAkBdBBAAAJEsQAQAAyRJEAABAsgQRAACQLEEEAAAkSxABAADJEkQAAECyBBEAAJAsQQQAXfbggw9GrVaLu+66K++hAHAIggiggLZv3x5XX311nHzyybFgwYJYtGhRnHPOOXHzzTfHr3/967yH19Gtt94an//85/MeBgDM2GF5DwCAVvfcc09cdtllMTQ0FJdffnmcfvrpMTIyEg899FB86EMfiscffzw2bdqU9zDbuvXWW+OYY46JK6+8Mu+hAMCMCCKAAtmxY0e89a1vjRUrVsQ3v/nN+N3f/d3Gfe95z3viySefjHvuuSfHEaYry7LYv39/HHHEEXkPBYAu8pY5gAL5xCc+EXv37o1//ud/bomhcatWrYo///M/b/z9pZdeio9//OOxcuXKGBoaiuHh4fjIRz4SL774Ystxw8PDcckll8RDDz0U69atiwULFsTJJ58cX/jCF6Y8xu7du+P9739/LF++PIaGhmLVqlVx4403xujo6LRjHx4ejscffzy+9a1vRa1Wi1qtFn/0R3/UuP+pp56Kyy67LI4++uj4rd/6rXj1q18947j7xje+Eeeee24sXrw4Fi5cGKecckp85CMfadnnxRdfjBtuuCFWrVoVQ0NDsXz58rjuuuumXIvbb789LrjggjjuuONiaGgoTj311PjMZz7T9ue55JJL4v77748zzzwzjjjiiPinf/qnxjX6i7/4ixgeHo6hoaFYtmxZXH755bFr166Wc4yOjsbf/M3fxLJly2LBggWxfv36ePLJJ2f0MwPQH7Usy7K8BwHAmGXLlsXQ0FBs3759RvtfeeWVcccdd8Sb3vSmOP/88+N73/tefOELX4hLL700vvKVrzT2Gx4ejgULFsTu3bvjT/7kT2Lp0qXxuc99Ln74wx/GY489FqeddlpERLzwwgtx9tlnxy9+8Yu4+uqr48QTT4zvfve78cUvfjGuueaa+NSnPtVxLHfffXe8733vi4ULF8ZHP/rRiIg4/vjj48ILL4ydO3fGmjVr4oUXXohrrrkmfud3fifuuOOOeOyxx+Kuu+6KN7zhDR3P+/jjj8erXvWqeOUrXxnveMc7YmhoKJ588sl45JFH4lvf+lZEjIXHa1/72njooYfi3e9+d6xevToee+yxuO2222Ljxo1x9913N863bt26OO2002LNmjVx2GGHxde+9rX4+te/Hrfccku85z3vablmhx9+ePzyl7+Mq6++OoaHh+OUU06JM888M84+++z4z//8z3jXu94Vr3rVq2LXrl2xZcuW2LRpU5xxxhnx4IMPxvnnnx+///u/H/V6Pd7+9rfHnj174hOf+EScfvrp8b3vfW9Gzy8AfZABUAh79uzJIiJ7/etfP6P9H3300Swisquuuqpl+wc/+MEsIrJvfvObjW0rVqzIIiLbtm1bY9uzzz6bDQ0NZR/4wAca2z7+8Y9nRx55ZPbTn/605ZzXX399NjAwkP385z+fdkynnXZadt55503Z/v73vz+LiOzb3/52Y9vzzz+fnXTSSdnw8HB24MCBjuf85Cc/mUVE9r//+78d9/niF7+Y1ev1lvNnWZbddtttWURk3/nOdxrbXnjhhSnHX3TRRdnJJ5/csm38mt13330t2//yL/8yi4js3/7t36acZ3R0NMuyLNu6dWsWEdnq1auzF198sXH/zTffnEVE9thjj3X8WQDoL2+ZAyiI5557LiIijjrqqBntf++990ZExLXXXtuy/QMf+EBExJS3o5166qnxmte8pvH3Y489Nk455ZR46qmnGts2b94cr3nNa2LJkiWxa9euxj8bNmyIAwcOxLZt22b/gx0c67p16+Lcc89tbFu4cGG8+93vjqeffjr+4z/+o+OxixcvjoiIr371qx3ftrd58+ZYvXp1vPzlL28Z9wUXXBAREVu3bm3s2/wZoD179sSuXbvivPPOi6eeeir27NnTct6TTjopLrroopZtX/7yl2PNmjVtV7VqtVrL39/5znfG4OBg4+/j17/5mgOQL0EEUBCLFi2KiIjnn39+Rvv/7Gc/i3q9HqtWrWrZfsIJJ8TixYvjZz/7Wcv2E088cco5lixZEv/3f//X+PsTTzwR9913Xxx77LEt/2zYsCEiIp599tlZ/UzNYz3llFOmbF+9enXj/k7e8pa3xDnnnBNXXXVVHH/88fHWt7417rzzzpY4euKJJ+Lxxx+fMu6XvexlU8b9ne98JzZs2BBHHnlkLF68OI499tjG55HaBdFk27dvj9NPP31GP/fka75kyZKIiJZrDkC+fMscQEEsWrQoli5dGj/+8Y9nddzkVYlOBgYG2m7Pmj5KOjo6GhdeeGFcd911bfcdD4x+OuKII2Lbtm2xdevWuOeee+K+++6LL33pS3HBBRfE17/+9RgYGIjR0dF4xSteETfddFPbcyxfvjwixmJm/fr18fKXvzxuuummWL58eQwODsa9994bn/zkJ6esQM33G+Vmcs0ByJcgAiiQSy65JDZt2hQPP/xwnH322dPuu2LFihgdHY0nnniisdISEbFz587YvXt3rFixYtaPv3Llyti7d29jRWi2OsXZihUr4r/+67+mbP/JT37SuH869Xo91q9fH+vXr4+bbrop/vZv/zY++tGPxtatW2PDhg2xcuXK+NGPfhTr16+fNhC/9rWvxYsvvhhbtmxpWb1pfkvdoaxcuXLW0QpAcXnLHECBXHfddXHkkUfGVVddFTt37pxy//bt2+Pmm2+OiIiLL744ImLKN7+Nr5Js3Lhx1o//5je/OR5++OG4//77p9y3e/fueOmll6Y9/sgjj4zdu3dP2X7xxRfHI488Eg8//HBj2759+2LTpk0xPDwcp556asdz/upXv5qy7YwzzoiIaHyl9pvf/Ob4xS9+EZ/97Gen7PvrX/869u3bFxETKzbNKzR79uyJ22+/fdqfq9kb3/jG+NGPftTyLX7jrPwAlI8VIoACWblyZfzLv/xLvOUtb4nVq1fH5ZdfHqeffnqMjIzEd7/73di8eXNceeWVERGxZs2auOKKK2LTpk2xe/fuOO+88+KRRx6JO+64Iy699NI4//zzZ/34H/rQh2LLli1xySWXxJVXXhlr166Nffv2Nb4e++mnn45jjjmm4/Fr166Nz3zmM/HXf/3XsWrVqjjuuOPiggsuiOuvvz7+9V//NV772tfGNddcE0cffXTccccdsWPHjvjyl78c9Xrn/3/ur/7qr2Lbtm2xcePGWLFiRTz77LNx6623xrJlyxpf0vCOd7wj7rzzzvizP/uz2Lp1a5xzzjlx4MCB+MlPfhJ33nln478l9Md//McxODgYr3vd6+Lqq6+OvXv3xmc/+9k47rjj4n/+539mfI3uuuuuuOyyy+Jd73pXrF27Nn71q1/Fli1b4rbbbos1a9bM7qIDkK98v+QOgHZ++tOfZn/6p3+aDQ8PZ4ODg9lRRx2VnXPOOdk//uM/Zvv372/s95vf/Cb72Mc+lp100knZ4Ycfni1fvjz78Ic/3LJPlo19hfTGjRunPM5555035Wuyn3/++ezDH/5wtmrVqmxwcDA75phjsj/8wz/M/v7v/z4bGRmZdtzPPPNMtnHjxuyoo47KIqLl3Nu3b8/e9KY3ZYsXL84WLFiQrVu3Lvv3f//3Q16LBx54IHv961+fLV26NBscHMyWLl2ave1tb5vy1eAjIyPZjTfemJ122mnZ0NBQtmTJkmzt2rXZxz72sWzPnj2N/bZs2ZK98pWvzBYsWJANDw9nN954Y/a5z30ui4hsx44dh7xmWZZlv/zlL7P3vve92e/93u9lg4OD2bJly7Irrrgi27VrV5ZlE1+7vXnz5pbjduzYkUVEdvvttx/y5wagP/yHWQEAgGT5DBEAAJAsQQQAACRLEAEAAMkSRAAAQLIEEQAAkCxBBAAAJEsQAQAAyRJEAABAsgQRAACQLEEEAAAk67C8BwBA71xYvyzvIVTCN0Y35z0EAHqklmVZlvcgAGhP0FSDoAIoLkEE0CfihtkQUQD9IYgA5kHkUATiCWDuBBFAB2KHKhFNAO0JIiBZggcmCCYgVYIIqCSxA90nmoAqEkRAaYkeKA6xBJSVIAIKTfRA+YkloMgEEVAIwgfSI5SAIhBEQF8JH+BQhBLQT4II6BnxA3SLSAJ6RRAB8yZ8gLwIJWC+BBEwK+IHKDqRBMyGIAI6Ej9AVYgkoBNBBESE+AHSI5KACEEEyRJAAK0EEqRJEEECxA/A3IgkqD5BBBUjfgB6SyRBtQgiKDkBBJAvgQTlJoigZAQQQLEJJCgXQQQFJ4AAyk0gQbEJIigYAQRQbQIJikUQQc4EEEDaBBLkSxBBnwkgAKYjkKC/BBH0gQgCYC7EEfSeIIIeEEAA9IJAgu4TRNAlIgiAfhJH0B2CCOZBBAFQBOII5k4QwSwIIADKQCDBzAkiOAQRBECZiSOYniCCNkQQAFUkjmAqQQQHiSAAUiKOYEw97wEAAADkxQoRSbMqBABWi0ibICI5IggAOhNHpEYQkQQRBACzJ45IgSCiskQQAHSPOKKqBBGVIYAAoH8EElUhiCg9IQQA+RFGlJ0gopREEAAUjziijAQRpSGCAKA8xBFlIYgoPCEEAOUljCg6QUQhiSAAqB5xRBEJIgpFCAFA9QkjikQQkTsRBADpEkfkTRCRGyEEAIwTRuRFENFXIggAOBRxRD8JIvpCCAEAsyWM6AdBRE8JIQBgvoQRvSSI6DoRBAD0ijii2wQRXSOEAIB+EUZ0iyBi3oQQAJAXYcR8CSLmTAgBAEUhjJgrQcSsCSEAoKiEEbMliJgxIQQAlIUwYqYEEdMSQQBA2YkjpiOIaEsIAQBVI4xoRxDRQggBAFUnjGgmiIgIIQQApEcYESGIkieEAIDUCaO0CaJECSEAgFbCKE31vAdA/4khAICpzJHSZIUoIX7JAQBmxmpROgRRAoQQAMDcCKPqE0QVJoQAALpDGFWXIKogIQQA0BvCqHoEUYUIIQCA/hBG1SGIKkAIAQD0nyiqBkFUYkIIAKAYxFF5CaISEkIAAMUkjMpHEJWIEAIAKAdhVB71vAfAzIghAIDyMHcrDytEBeeXCQCg3KwWFZsgKighBABQLcKomLxlroDEEABA9ZjjFZMVogLxSwIAkAarRcVhhaggxBAAQDrM/YrDClHO/DIAAKTNalG+rBDlSAwBAGBOmC8rRDnwogcAoB2rRf1nhajPxBAAAJ2YK/afFaI+8eIGAGA2rBb1hxWiPhBDAADMljlkf1gh6iEvYgAAusFqUe9YIeoRMQQAQLeYW/aOFaIu82IFAKCXrBZ1lxWiLhJDAAD0mjlndwmiLvHCBACgX8w9u8db5ubJixEAgDx5C938WCGaBzEEAEDezEnnRxDNkRceAABFYW46d94yN0tebAAAFJm30M2OFaJZEEMAABSdOevsCKIZ8sICAKAszF1nThDNgBcUAABlYw47Mz5DNA0vIgAAqsDnijqzQtSBGAIAoCrMbTsTRG14wQAAUDXmuO0Jokm8UAAAqCpz3akEEQAAkCxB1EQxAwBQdea8rQTRQV4YAACkwtx3QvJfu+3FAABAylL/Su6kV4jEEAAAqUt9TpxsEKX+xAMAwLiU58ZJBlHKTzgAALST6hw5uSBK9YkGAIBDSXGunFQQpfgEAwDAbKQ2Z04miFJ7YgEAYK5SmjsnEUQpPaEAANANqcyhKx9EqTyRAADQbSnMpSsdRCk8gQAA0EtVn1NXNoiq/sQBAEC/VHluXckgqvITBgAAeajqHLtyQVTVJwoAAPJWxbl2pYKoik8QAAAUSdXm3JUJoqo9MQAAUFRVmntXIoiq9IQAAEAZVGUOXvogqsoTAQAAZVOFuXipg6gKTwAAAJRZ2efkpQ4iAACA+ahlWZblPYjZKnuFAgBAFX1jdHPeQ5i10q0QiSEAACimMs7VSxVEZbzAAACQkrLN2UsTRGW7sAAAkKoyzd1LEURluqAAAEB55vClCCIAAIBeKHwQlaUsAQCAVmWYyxc6iMpwAQEAgM6KPqcvbBAV/cIBAAAzU+S5fWGDCAAAoNcKGURFLkgAAGD2ijrHL1wQFfVCAQAA81PEuX6hgqiIFwgAAOieos35CxVEAAAA/VSYICpaKQIAAL1RpLl/IYKoSBcEAADovaI0QCGCCAAAIA+5B1FRyhAAAOivIrRArkFUhAsAAADkJ+8myH2FCAAAIC+5BVHeJQgAABRDnm2QSxCJIQAAoFlejeAtcwAAQLL6HkRWhwAAgHbyaAUrRAAAQLL6GkRWhwAAgOn0uxmsEAEAAMnqWxBZHQIAAGain+1ghQgAAEhWX4LI6hAAADAb/WqIngeRGAIAAOaiHy3hLXMAAECyBBEAAJCsngaRt8sBAADz0eumsEIEAAAkq2dBZHUIAADohl62hRUiAAAgWT0JIqtDAABAN/WqMawQAQAAyep6EFkdAgAAeqEXrWGFCAAASJYgAgAAktXVIPJ2OQAAoJe63RxWiAAAgGQJIgAAIFldCyJvlwMAAPqhm+1hhQgAAEiWIAIAAJIliAAAgGR1JYh8fggAAOinbjWIFSIAACBZgggAAEjWvIPI2+UAAIA8dKNFrBABAADJEkQAAECyBBEAAJCseQWRzw8BAAB5mm+TWCECAACSJYgAAIBkCSIAACBZcw4inx8CAACKYD5tYoUIAABIliACAACSJYgAAIBkCSIAACBZgggAAEjWnILIN8wBAABFMtdGsUIEAAAkSxABAADJEkQAAECyBBEAAJAsQQQAACRLEAEAAMkSRAAAQLIEEQAAkCxBBAAAJEsQAQAAyRJEAABAsgQRAACQLEEEAAAkSxABAADJEkQAAECyBBEAAJAsQQQAACRLEAEAAMkSRAAAQLIEEQAAkCxBBAAAJEsQAQAAyRJEAABAsgQRAACQLEEEAAAkSxABAADJEkQAAECy5hRE3xjd3O1xAAAAzNlcG8UKEQAAkCxBBAAAJEsQAQAAyRJEAABAsgQRAACQrDkHkW+aAwAAimA+bWKFCAAASJYgAgAAkiWIAACAZM0riHyOCAAAyNN8m8QKEQAAkCxBBAAAJEsQAQAAyZp3EPkcEQAAkIdutIgVIgAAIFmCCAAASFZXgsjb5gAAgH7qVoNYIQIAAJIliAAAgGQJIgAAIFldCyKfIwIAAPqhm+1hhQgAAEiWIAIAAJLV1SDytjkAAKCXut0cVogAAIBkCSIAACBZXQ8ib5sDAAB6oRetYYUIAABIVk+CyCoRAADQTb1qDCtEAABAsnoWRFaJAACAbuhlW1ghAgAAktXTILJKBAAAzEevm8IKEQAAkCxBBAAAJKvnQeRtcwAAwFz0oyX6skIkigAAgNnoV0N4yxwAAJCsvgWRVSIAAGAm+tkOVogAAIBk9TWIrBIBAADT6XczWCECAACS1fcgskoEAAC0k0crWCECAACSlUsQWSUCAACa5dUIua0QiSIAACAi3zbwljkAACBZuQaRVSIAAEhb3k2Q+wpR3hcAAADIRxFaIPcgAgAAyEshgqgIZQgAAPRPURqgEEEUUZwLAgAA9FaR5v6FCSIAAIB+K1QQFakUAQCA7ivanL9QQRRRvAsEAAB0RxHn+oULoohiXigAAGDuijrHL2QQAQAA9ENhg6ioBQkAAMxOkef2hQ2iiGJfOAAA4NCKPqcvdBBFFP8CAgAA7ZVhLl/4IAIAAOiVUgRRGcoSAACYUJY5fCmCKKI8FxQAAFJXprl7aYIoolwXFgAAUlS2OXupgiiifBcYAABSUca5ei3LsizvQczVhfXL8h4CAAAkr4whNK50K0QAAADdUuogKnOJAgBAFZR9Tl7qIIoo/xMAAABlVYW5eOmDKKIaTwQAAJRJVebglQiiiOo8IQAAUHRVmntXJogiqvXEAABAEVVtzl2pIIqo3hMEAABFUcW5duWCKKKaTxQAAOSpqnPsSgZRRHWfMAAA6Lcqz60rG0QR1X7iAACgH6o+p650EEVU/wkEAIBeSWEuXfkgikjjiQQAgG5KZQ6dRBBFpPOEAgDAfKU0d04miCLSemIBAGAuUpszJxVEEek9wQAAMFMpzpWTC6KINJ9oAACYTqpz5CSDKCLdJxwAACZLeW6cbBBFpP3EAwBAhDlxLcuyLO9BFMGF9cvyHgIAAPRN6iE0LukVomZeEAAApMLcd4IgAgAAkiWImihlAACqzpy3lSCaxAsEAICqMtedShC14YUCAEDVmOO2J4g68IIBAKAqzG0787XbM+AruQEAKCMhdGhWiGbACwkAgLIxh50ZQTRDXlAAAJSFuevMCaJZ8MICAKDozFlnx2eI5sjnigAAKBIhNDdWiObICw4AgKIwN507QTQPXngAAOTNnHR+vGWuS7yFDgCAfhJC3WGFqEu8IAEA6Bdzz+4RRF3khQkAQK+Zc3aXt8z1iLfQAQDQTUKoN6wQ9YgXLAAA3WJu2TtWiPrAahEAAHMhhHrPClEfeCEDADBb5pD9YYWoz6wWAQAwHSHUX1aI+swLHACATswV+88KUY6sFgEAECGE8mSFKEde+AAAmBPmywpRQVgtAgBIixAqBitEBeEXAgAgHeZ+xWGFqICsFgEAVJMQKh4rRAXkFwUAoHrM8YrJClHBWS0CACg3IVRsgqgkhBEAQLkIoXLwlrmS8AsFAFAe5m7lYYWohKwWAQAUkxAqH0FUYsIIAKAYhFB5CaIKEEYAAP0ngqpBEFWIMAIA6A8xVB2CqIKEEQBAbwih6hFEFSaMAAC6QwhVlyBKgDACAJgbIVR9gighwggAYGaEUDoEUYKEEQBAe0IoPfW8B0D/+UUHAJjKHClNVogSZ7UIAEidEEqbICIihBEAkB4hRIQgYhJhBABUnRCimSCiLWEEAFSNEKIdQcS0hBEAUHZCiOkIImZMHAEAZSGCmClBxKwJIwCgqIQQsyWImDNhBAAUhRBirgQR8yaMAIC8CCHmSxDRNcIIAOgXIUS3CCK6ThgBAL0ihOg2QURPiSMAYL5EEL0kiOgLYQQAzJYQoh8EEX0ljACAQxFC9JMgIjfiCAAYJ4LIiyAid8IIANIlhMibIKJQxBEAVJ8IokgEEYUkjACgeoQQRSSIKDxxBADlJYIoOkFEaQgjACgPIURZCCJKSRwBQPGIIMpIEFF64ggA8iOCKDtBRGUIIwDoHyFEVQgiKksgAUD3CCCqShCRBHEEALMngkiBICI54ggAOhNBpEYQkTRxBAAiiLQJIjhIHAGQEhEEY+p5DwAAACAvVoigDatFAFSRVSGYShDBIYgjAMpMBMH0BBHMgjgCoAxEEMycIIJ5EEgAFIEAgrkTRNAl4giAfhJB0B2CCHpAHJGy+//fo3HR0jNa/j4fk8/V/HdIjQiC7hNE0AcCiaqZb+T0iliiagQQ9J4ggj4TR5RNUeNnpkQSZSOCoL8EEeRMIFE0ZQ+gQxFIFI0AgnwJIigYgUSeqh5D40QReRJAUCyCCApOINFLqQTQoQgkekkAQbEJIigZgcR8iaDpiSPmSwBBuQgiKDmBxEyIoLkRR8yEAIJyE0RQMQKJycTQ/IgiJhNAUC2CCBIgktIkhLpLGKVJ/ED1HZb3AIDea/cvdJEE0Er8QJqsEAERIZCqxupQb1glqhYBBEQIImAaIqmcxFBviaJyEj9AJ4IImBWRVGxiqD9EUbGJH2A2fIYImBWfRwKKRPwA8yWIgHnrNCERSkC3CB+gVwQR0DNWk4C5ED9AP/kMEVAIQql7fI6ot3x+qHuED1AEVoiAQvC2O6gu4QMUmRUioLTEUmdWiXrD6lBnogcoK0EEVJJYGiOMukMIjRE9QBUJIiBZqUSTKJqfVGJI7ACpEkQAHVQxmMTRzFQxggQPQHuCCGAeyhpNwqi9soaQ2AGYO0EE0CdFj6fUIqno8SNyAPpDEAEUWJ4RVbVAyjOAxA1AcQkigArrRVAVPZR6ET6CBqC6BBEAAJCset4DAAAAyIsgAgAAkiWIAACAZAkiAAAgWYIIAABIliACAACSJYgAAIBkCSIAACBZgggAAEiWIAIAAJIliAAAgGQJIgAAIFmCCAAASJYgAgAAkiWIAACAZAkiAAAgWYIIAABIliACAACSJYgAAIBkCSIAACBZgggAAEiWIAIAAJIliAAAgGQJIgAAIFmCCAAASNZheQ8Aqmj//v0xMjKS9zAAqJjBwcFYsGBB3sOAShFE0GX79++P3z5iSYzE/ryHAkDFnHDCCbFjxw5RBF0kiKDLRkZGYiT2x7lxcRxWG4pavTZ2R63edPvgn/Va43atXm/aXp/Yb/z+Wn1s/+bja7VJ+0bT9sn7Tj6+81iyWm3iDbXNjzXt/Qe3H7zdsu/kbfXm+5vOc/B28/mzlu0T+2YxfjsaP9fE/ZP2bdre2K/58ce319sf39ByfIfbba5Fy+O33bfN7ehw/6SxdDy+07ZDPP64jtvanKvdtYhaNuOxRG3sp2r8bFP2zdo/ZvP2xuNPbKtNd3xkTb8CE49c63D8+PaWczYdX5t0fNOvS9Rbjp+4v960rX7wp28+z/i+9Um3IyLq0bqt3ub2+Lk63T/+mGPbRiceKybfPxoDTcdM7Dt2/oHIotZ0/MS+Tduabx/cd/xxBmqjjXMOHHy88fM2HqvNuQZqo40xDjTtN/5rPBDN5x0/Jmuca+z48eMmzjPQ9PMPNI1l/HlpnCsmruX4/WPbYuJajY+lFjFw8BmZ2FaLemPbxO2B2vi2etO2sdvPPT8aK9Y+HSMjI4IIukgQQY8cFofHYbXDo9YUJM23x/6cmDXVak1B1Hx/ven+QwZRbcpxHYOoNun+eQdRbSICWoJnIgy6GkTNE+MeB1H74zvcbhxfa7rddPwhguRQEdHzIGpz/7jJQTTddZlTEDU9Vt+DqPl2m+MnB9HkYyZeLhMT5zkFUZtj5hpErcEz8yBq3j72Z6cgao6A2QdRvW0QZR1uzzyIBhrjqkX94AUdj5GxIBq/XWuKlKxpW9Z0rmiMZeLxJ7ZNF0QDswiigRkFkY99Q6/47QIAAJIliAAAgGQJIgAAIFmCCAAASJYgAgAAkiWIAACAZAkiAAAgWYIIAABIliACAACSJYgAAIBkCSIAACBZgggAAEjWYXkPAKrqpfhNRFaPWlY7uKX59sE/s1rjdi2rN22vT+w3evD+WtP9teY/6xO3G6evtdl38vGT7s9qjdtZrRaRxaTjD3X/we0REaPRuu/4jzy+rd58f9N5GkOpNX6srGX7xL7jDx+1iKhPPtekfdte9qZj2lz2lnFFu+M73J58LWqt29vv2+Z2dLh/0lg6Ht9p2yEef1zHbW3O1e5aRC2b8ViiNvZTNX62Kftm7R+zeXvj8Se21aY7PrKmX4GJR651OH58e8s5m46vTTq+6dclspbjJ+7PmrZlB3/68fOMNh1Tn3Q7IqIerdvqbW43/penw/31aN42OvFYMfn+0RhoOmZi37HzD0QWtabjJ/Zt2tZ8++C+448zUBttnHPg4OONn7fxWG3ONVAbbYxxoGm/8V/jgWg+7/gxWeNcY8ePHzdxnoGmn3+gaSzjz1XjXDFxLcfvH9sWE9dqfCy1iIGDz8jEtlrUG9smbg/UJvab2DZ2vueeHw2g+wQRdFmWZbFw4cJ4aO+9YzO8A3mPCICqWLhwYWRZdugdgRkTRNBltVot9u7dG//93/8dixYtyns4AFTEc889F8uXL49a88o1MG+CCHpk0aJFgggAoOB8qQIAAJAsQQQAACRLEEGXDQ0NxQ033BBDQ0N5DwWACvHvF+iNWuarSgAAgERZIQIAAJIliAAAgGQJIgAAIFmCCAAASJYgAgAAkiWIoIs+/elPx/DwcCxYsCDOOuuseOSRR/IeEgAlt23btnjd614XS5cujVqtFnfffXfeQ4JKEUTQJV/60pfi2muvjRtuuCF+8IMfxJo1a+Kiiy6KZ599Nu+hAVBi+/btizVr1sSnP/3pvIcCleS/QwRdctZZZ8Uf/MEfxC233BIREaOjo7F8+fJ43/veF9dff33OowOgCmq1WnzlK1+JSy+9NO+hQGVYIYIuGBkZie9///uxYcOGxrZ6vR4bNmyIhx9+OMeRAQAwHUEEXbBr1644cOBAHH/88S3bjz/++HjmmWdyGhUAAIciiAAAgGQJIuiCY445JgYGBmLnzp0t23fu3BknnHBCTqMCAOBQBBF0weDgYKxduzYeeOCBxrbR0dF44IEH4uyzz85xZAAATOewvAcAVXHttdfGFVdcEWeeeWasW7cuPvWpT8W+ffvine98Z95DA6DE9u7dG08++WTj7zt27IhHH300jj766DjxxBNzHBlUg6/dhi665ZZb4u/+7u/imWeeiTPOOCP+4R/+Ic4666y8hwVAiT344INx/vnnT9l+xRVXxOc///n+DwgqRhABAADJ8hkiAAAgWYIIAABIliACAACSJYgAAIBkCSIAACBZgggAAEiWIAIAAJIliAAAgGQJIgAAIFmCCAAASJYgAgAAkvX/AWxtjpnl6gl2AAAAAElFTkSuQmCC" + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "## Plot the cone using healpy for demonstration\n", + "\n", + "NSIDE = 256\n", + "NPIX = hp.nside2npix(NSIDE)\n", + "m = np.zeros(NPIX)\n", + "center_vec = hp.ang2vec(ra, dec, lonlat=True)\n", + "radius_radians = np.radians(radius)\n", + "cone_pixels = hp.query_disc(NSIDE, center_vec, radius_radians)\n", + "m[cone_pixels] = 1\n", + "hp.mollview(m, title=\"Cone to search\")" + ], + "metadata": { + "collapsed": false + } + }, + { + "cell_type": "code", + "execution_count": 35, + "outputs": [ + { + "data": { + "text/plain": "
", + "image/png": "iVBORw0KGgoAAAANSUhEUgAAA0QAAAIECAYAAAA5Nu72AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/P9b71AAAACXBIWXMAAA9hAAAPYQGoP6dpAAA1uUlEQVR4nO3dd5hcdb348c/sJrsJaYBpgJiQoFKCwAWkRA0XcgktBFQQkJCAIEqogvx4VAyxwKUoUqRZwCuKFCECgihVafcqTUHpBAUlEAihmBDInt8fYWZntmXL7MyZOa/X8/AwOXPOzJmZs/B953vmbC5JkiQAAAAyqKHaOwAAAFAtgggAAMgsQQQAAGSWIAIAADJLEAEAAJkliAAAgMwSRAAAQGYJIgAAILMEEQAAkFmCCOiRBQsWRC6Xi0svvbTau9KpSuzj+PHjY/bs2b3aNk3vYS6Xi5NPPrnau0EfXXrppZHL5WLBggWFZdtvv31sv/32ZX2e/LF75plnlvVxq8GxD+QJIki5p59+Og477LCYMGFCDBo0KIYPHx6TJ0+Os88+O5YuXdrjxzv//PNTMRAnne655544+eST47XXXqv2rkBN+va3vx177LFHjBkzRnRBjRhQ7R0AOvfrX/869t5772hubo4DDzwwJk2aFMuXL4+77rorvvzlL8ejjz4aF198cY8e8/zzz4+RI0f2enajFowbNy6WLl0aAwcOrPaupN7SpUtjwIDW/xXcc889MW/evJg9e3asvvrq1dsxqFFf+9rXYuzYsbH55pvHzTffXO3dAbpBEEFKPfvss7HvvvvGuHHj4rbbbou11lqrcN+cOXPiqaeeil//+tdV3MP0yuVyMWjQoGrvRk3wPlGv3nrrrRgyZEhFnmvZsmXR1NQUDQ0N8eyzz8b48eNj0aJFMWrUqIo8P9A3TpmDlDr99NPjzTffjB/96EclMZS3/vrrx9FHH1348yWXXBI77LBDjB49Opqbm2OjjTaKCy64oGSb8ePHx6OPPhp33nln5HK5yOVyhe8YvPrqq3H88cfHJptsEkOHDo3hw4fHLrvsEg8//HC39ve2226Lj3/84zFkyJBYffXVY8aMGfG3v/2t3Xp33HFHbLnlljFo0KCYOHFiXHTRRXHyySdHLpdb5XNsv/32MWnSpLj//vtju+22i8GDB8d6660XF154Ycl6bb+j89JLL8WoUaNi++23jyRJCus99dRTMWTIkPjMZz5TWPb222/H3LlzY/3114/m5uZYd91144QTToi33367W+9DW6+99lrMnj07RowYEauvvnrMmjWr09PRHnvssfj0pz8da665ZgwaNCi23HLLuO6660rWyX9X5O67744vfelLMWrUqBgyZEjstdde8fLLL5es+6c//SmmTZsWI0eOLLxXBx98cMk6xaf0nHzyyfHlL385IiLWW2+9wjGyYMGCmDJlSmy66aYd7veHP/zhmDZtWpfvw/jx42P33XcvfP6DBw+OTTbZJO64446IiLjmmmtik002iUGDBsUWW2wRDz74YMn2f/7zn2P27NmFU0fHjh0bBx98cLzyyisl6+WPpcceeyz22WefGD58eLzvfe+Lo48+OpYtW9blPvbUiy++GAcddFC8//3vj+bm5lhrrbVixowZJd/jqdTrLpfuHDNtJUkSn//856OpqSmuueaaPh8rbV111VWxxRZbxODBg2PkyJFxwAEHxAsvvFCyzuzZs2Po0KHx9NNPx6677hrDhg2Lz372sxGx8mf62GOPjVGjRsWwYcNijz32iOeff77D53rhhRfi4IMPjjFjxkRzc3NsvPHG8eMf/7hknTvuuCNyuVz84he/iK997WuxzjrrxGqrrRavv/56RKz8zIHaYoYIUur666+PCRMmxHbbbdet9S+44ILYeOONY4899ogBAwbE9ddfH4cffni0tLTEnDlzIiLie9/7Xhx55JExdOjQ+OpXvxoREWPGjImIiGeeeSbmz58fe++9d6y33nqxcOHCuOiii2LKlCnx17/+NdZee+1On/uWW26JXXbZJSZMmBAnn3xyLF26NM4999yYPHlyPPDAA4UBwoMPPhg777xzrLXWWjFv3rxYsWJFfOMb3+jR36IuXrw4dt1119hnn31iv/32iyuvvDK++MUvRlNTU6cDt9GjR8cFF1wQe++9d5x77rlx1FFHRUtLS8yePTuGDRsW559/fkREtLS0xB577BF33XVXfP7zn48NN9ww/vKXv8RZZ50VTzzxRMyfP7/b+xmxcqA4Y8aMuOuuu+ILX/hCbLjhhnHttdfGrFmz2q376KOPxuTJk2OdddaJE088MYYMGRJXXnll7LnnnvHLX/4y9tprr5L1jzzyyFhjjTVi7ty5sWDBgvje974XRxxxRFxxxRURsTICd9pppxg1alSceOKJsfrqq8eCBQvimmuu6XR/P/nJT8YTTzwRl19+eZx11lkxcuTIiIgYNWpUzJw5Mw499NB45JFHYtKkSYVt/vjHP8YTTzwRX/va11b5fjz11FOx//77x2GHHRYHHHBAnHnmmTF9+vS48MIL4ytf+UocfvjhERFx6qmnxj777BOPP/54NDSs/Hu73/3ud/HMM8/EQQcdFGPHji2cLvroo4/Gfffd1y6o99lnnxg/fnyceuqpcd9998U555wTixcvjv/5n/9Z5X5216c+9al49NFH48gjj4zx48fHSy+9FL/73e/i73//e8mguJKvuy96c8ysWLEiDj744Ljiiivi2muvjd122y1effXVPh8reZdeemkcdNBBsdVWW8Wpp54aCxcujLPPPjvuvvvuePDBB0tO63z33Xdj2rRp8bGPfSzOPPPMWG211SIi4pBDDonLLrss9t9//9huu+3itttui912263dcy1cuDC22WabyOVyccQRR8SoUaPipptuis997nPx+uuvxzHHHFOy/je/+c1oamqK448/Pt5+++1oamrq9usCUiYBUmfJkiVJRCQzZszo9jb//ve/2y2bNm1aMmHChJJlG2+8cTJlypR26y5btixZsWJFybJnn302aW5uTr7xjW+ULIuI5JJLLiks22yzzZLRo0cnr7zySmHZww8/nDQ0NCQHHnhgYdn06dOT1VZbLXnhhRcKy5588slkwIABSXf+czRlypQkIpLvfOc7hWVvv/124fmXL1/e6T4mSZLst99+yWqrrZY88cQTyRlnnJFERDJ//vzC/T/96U+ThoaG5A9/+EPJdhdeeGESEcndd99dWDZu3Lhk1qxZXe7v/Pnzk4hITj/99MKyd999N/n4xz/ebv923HHHZJNNNkmWLVtWWNbS0pJst912yQc/+MHCsksuuSSJiGTq1KlJS0tLYfmxxx6bNDY2Jq+99lqSJEly7bXXJhGR/PGPf+xyHyMimTt3buHP+ffl2WefLVnvtddeSwYNGpT8v//3/0qWH3XUUcmQIUOSN998s8vnGTduXBIRyT333FNYdvPNNycRkQwePDh57rnnCssvuuiiJCKS22+/vbCso+P78ssvTyIi+f3vf19YNnfu3CQikj322KNk3cMPPzyJiOThhx/ucj+7a/HixUlEJGeccUaX61XqdeePi+LPbcqUKR3+rHemO8dM/mfrjDPOSN55553kM5/5TDJ48ODk5ptvLqzT12Mlb/ny5cno0aOTSZMmJUuXLi0sv+GGG5KISL7+9a8Xls2aNSuJiOTEE08seYyHHnooiYjk8MMPL1m+//77tzv2P/e5zyVrrbVWsmjRopJ1991332TEiBGFz+L2229PIiKZMGFCh59P3ssvv9zuOYB0csocpFD+1Ithw4Z1e5vBgwcXbi9ZsiQWLVoUU6ZMiWeeeSaWLFmyyu2bm5sLfyu9YsWKeOWVV2Lo0KHx4Q9/OB544IFOt/vXv/4VDz30UMyePTvWXHPNwvKPfOQj8V//9V9x4403Fh7zlltuiT333LNktmn99dePXXbZpduvc8CAAXHYYYcV/tzU1BSHHXZYvPTSS3H//fd3ue15550XI0aMiE9/+tNx0kknxcyZM2PGjBmF+6+66qrYcMMNY4MNNohFixYV/tlhhx0iIuL222/v9n5GRNx4440xYMCA+OIXv1hY1tjYGEceeWTJeq+++mrcdtttsc8++8Qbb7xReN5XXnklpk2bFk8++WS7U4Q+//nPl8wOfPzjH48VK1bEc889FxFR+JvzG264Id55550e7XdHRowYETNmzIjLL7+8cNrhihUr4oorrog999yzW9/V2GijjWLbbbct/HnrrbeOiIgddtghPvCBD7Rb/swzzxSWFR/fy5Yti0WLFsU222wTEdHh8ZmfFc3Lv+f547GvBg8eHE1NTXHHHXfE4sWLu1y3kq+7L3pyzCxfvjz23nvvuOGGG+LGG2+MnXbaqXBfOY6ViJWn77300ktx+OGHl3zXbbfddosNNtigw+9QFv+sRbR+3kcddVTJ8razPUmSxC9/+cuYPn16JElS8vM/bdq0WLJkSbv3e9asWSWfD1C7BBGk0PDhwyMi4o033uj2NnfffXdMnTq18B2eUaNGxVe+8pWIiG4FUUtLS5x11lnxwQ9+MJqbm2PkyJExatSo+POf/9zl9vkB+Ic//OF292244YaxaNGieOutt+Kll16KpUuXxvrrr99uvY6WdWbttdduN6D60Ic+FBFR8t2Njqy55ppxzjnnxJ///OcYMWJEnHPOOSX3P/nkk/Hoo4/GqFGjSv7JP/5LL73U7f2MWPnerLXWWjF06NCS5W3fq6eeeiqSJImTTjqp3XPPnTu3w+cuHkhHRKyxxhoREYXB+ZQpU+JTn/pUzJs3L0aOHBkzZsyISy65pNffhYqIOPDAA+Pvf/97/OEPf4iIladKLly4MGbOnNmt7dvu84gRIyIiYt111+1weXFovPrqq3H00UfHmDFjYvDgwTFq1KhYb731IqLj4/uDH/xgyZ8nTpwYDQ0NXR4jb775Zrz44ouFf9p+J6tYc3NznHbaaXHTTTfFmDFj4hOf+EScfvrp8eKLL1b1dfdFT46ZU089NebPnx9XX311h7/rqK/HSkTX/23ZYIMNCvfnDRgwIN7//ve3e4yGhoaYOHFiyfK2j/nyyy/Ha6+9FhdffHG7n8GDDjooItr/DOY/B6D2+Q4RpNDw4cNj7bXXjkceeaRb6z/99NOx4447xgYbbBDf/e53Y911142mpqa48cYb46yzzoqWlpZVPsYpp5wSJ510Uhx88MHxzW9+M9Zcc81oaGiIY445plvb14r8ZXAXL14czz//fMl3EFpaWmKTTTaJ7373ux1u23YAWy759/f444/v9AvnbaOxsbGxw/XyfyOfy+Xi6quvjvvuuy+uv/76uPnmm+Pggw+O73znO3Hfffe1i7TumDZtWowZMyYuu+yy+MQnPhGXXXZZjB07NqZOndqt7Tvb51W9loiV3wm655574stf/nJsttlmMXTo0GhpaYmdd965W8dnd75rc+aZZ8a8efMKfx43blyXAXXMMcfE9OnTY/78+XHzzTfHSSedFKeeemrcdtttsfnmm6/y9VXidfdET46ZadOmxW9+85s4/fTTY/vtt293tcK+Hiu9UTzL3VP59/KAAw7o8Dt+EStnvYuZHYL6IYggpXbfffe4+OKL49577y053aYj119/fbz99ttx3XXXlfxtdEeneHU2MLz66qvjP//zP+NHP/pRyfLXXnut8OX6jowbNy4iIh5//PF29z322GMxcuTIGDJkSAwaNCgGDRoUTz31VLv1OlrWmX/+85/tLqf7xBNPRMSqr+70m9/8Jn74wx/GCSecED/72c9i1qxZ8b//+7+F38MzceLEePjhh2PHHXcsy5fVx40bF7feemu8+eabJYPJtu/VhAkTIiJi4MCBZR8wbrPNNrHNNtvEt7/97fj5z38en/3sZ+MXv/hFHHLIIR2u39XrbmxsjP333z8uvfTSOO2002L+/Plx6KGHdjqwL5fFixfHrbfeGvPmzYuvf/3rheVPPvlkp9s8+eSTJX+D/9RTT0VLS0uXx8iBBx4YH/vYxwp/7s6Ad+LEiXHcccfFcccdF08++WRsttlm8Z3vfCcuu+yyVW67Kr153eXQnWNmm222iS984Qux++67x9577x3XXnttye+zKsexUvzflvxpq3mPP/544f5VPUZLS0s8/fTTJbNCbX8G81egW7FiRb9GG5BOTpmDlDrhhBNiyJAhccghh8TChQvb3f/000/H2WefHRGtf9Nc/DfLS5YsiUsuuaTddkOGDOnwss+NjY0l20es/E5N2++utLXWWmvFZpttFj/5yU9KHveRRx6J3/72t7HrrrsWHn/q1Kkxf/78+Oc//1lY76mnnoqbbrqpy+co9u6778ZFF11U+PPy5cvjoosuilGjRsUWW2zR6XavvfZaHHLIIfHRj340TjnllPjhD38YDzzwQJxyyimFdfbZZ5944YUX4gc/+EG77ZcuXRpvvfVWt/czImLXXXeNd999t+Ty5ytWrIhzzz23ZL3Ro0fH9ttvHxdddFH861//avc4XZ261ZnFixe3+zw322yziIguT5vLh2ZnlwafOXNmLF68OA477LB4880344ADDujxvvVUR8d3xMqrJnbm+9//fsmf8+95V99XmzBhQkydOrXwz+TJkztd99///ne7y3hPnDgxhg0b1qfTEov15nX3RU+PmalTp8YvfvGL+M1vfhMzZ85sN2PV12Nlyy23jNGjR8eFF15Y8vw33XRT/O1vf+vwSnFt5T/vtqfHtn0PGxsb41Of+lT88pe/7HBmvjc/g0DtMEMEKTVx4sT4+c9/Hp/5zGdiww03jAMPPDAmTZoUy5cvj3vuuSeuuuqqmD17dkRE7LTTTtHU1BTTp08vDD5+8IMfxOjRo9sNsLfYYou44IIL4lvf+lasv/76MXr06Nhhhx1i9913j2984xtx0EEHxXbbbRd/+ctf4mc/+1lh9qIrZ5xxRuyyyy6x7bbbxuc+97nCZbdHjBhR+B03ESt/R8xvf/vbmDx5cnzxi1+MFStWxHnnnReTJk2Khx56qFvvy9prrx2nnXZaLFiwID70oQ/FFVdcEQ899FBcfPHFMXDgwE63O/roo+OVV16JW265JRobG2PnnXeOQw45JL71rW/FjBkzYtNNN42ZM2fGlVdeGV/4whfi9ttvj8mTJ8eKFSviscceiyuvvDJuvvnm2HLLLbu1nxER06dPj8mTJ8eJJ54YCxYsiI022iiuueaaDr/78f3vfz8+9rGPxSabbBKHHnpoTJgwIRYuXBj33ntvPP/8893+fVB5P/nJT+L888+PvfbaKyZOnBhvvPFG/OAHP4jhw4cXIrUj+aj86le/Gvvuu28MHDgwpk+fXgilzTffPCZNmlS4AMV//Md/9Gi/emP48OGF7+i88847sc4668Rvf/vbePbZZzvd5tlnn4099tgjdt5557j33nsLl13u7Pfj9NQTTzwRO+64Y+yzzz6x0UYbxYABA+Laa6+NhQsXxr777luW5+jN6+6L3hwze+65Z1xyySVx4IEHxvDhw0v+sqKvx8rAgQPjtNNOi4MOOiimTJkS++23X+Gy2+PHj49jjz12lY+x2WabxX777Rfnn39+LFmyJLbbbru49dZbO5yV/u///u+4/fbbY+utt45DDz00Ntpoo3j11VfjgQceiFtuuSVeffXVbu33T3/603juuefi3//+d0RE/P73v49vfetbEbEyErszswVUWDUubQd03xNPPJEceuihyfjx45OmpqZk2LBhyeTJk5Nzzz235BLN1113XfKRj3wkGTRoUDJ+/PjktNNOS3784x+3uxTviy++mOy2227JsGHDkogoXJZ32bJlyXHHHZestdZayeDBg5PJkycn9957b7tL93Z2SetbbrklmTx5cjJ48OBk+PDhyfTp05O//vWv7V7Prbfemmy++eZJU1NTMnHixOSHP/xhctxxxyWDBg1a5XsxZcqUZOONN07+9Kc/Jdtuu20yaNCgZNy4ccl5551Xsl7bffzVr37V7nLdSZIkr7/+ejJu3Lhk0003LVyye/ny5clpp52WbLzxxklzc3OyxhprJFtssUUyb968ZMmSJYVtu3PZ7SRJkldeeSWZOXNmMnz48GTEiBHJzJkzkwcffLDD9/Dpp59ODjzwwGTs2LHJwIEDk3XWWSfZfffdk6uvvrqwTv7yym0vjZy/FHD+ks0PPPBAst9++yUf+MAHkubm5mT06NHJ7rvvnvzpT38q2S46uCzwN7/5zWSdddZJGhoaOrwE9+mnn55ERHLKKaes8vXnjRs3Ltltt93aLY+IZM6cOSXLii/tnPf8888ne+21V7L66qsnI0aMSPbee+/kn//8Z7v9z192+69//Wvy6U9/Ohk2bFiyxhprJEcccUTJpZv7atGiRcmcOXOSDTbYIBkyZEgyYsSIZOutt06uvPLKqrzuclx2uzvHTEf7mCRJcv755ycRkRx//PEly3tzrLR1xRVXJJtvvnnS3NycrLnmmslnP/vZ5Pnnny9ZZ9asWcmQIUM63H7p0qXJUUcdlbzvfe9LhgwZkkyfPj35xz/+0eGxv3DhwmTOnDnJuuuumwwcODAZO3ZssuOOOyYXX3xxYZ38z9pVV13V4fPlfz1AR/8UX1IdSI9ckrSZHweosD333DMeffTRVX43Yvvtt49FixZ1+2IT9I+zzz47jj322FiwYEG7K6hV28knnxzz5s2Ll19+ucvvvlEZaT5WAPJ8hwioqKVLl5b8+cknn4wbb7yxw0v3kj5JksSPfvSjmDJligEuXXKsALXCd4iAipowYULMnj07JkyYEM8991xccMEF0dTUFCeccEK1d40uvPXWW3HdddfF7bffHn/5y1/iV7/6VbV3iR54+eWXY8WKFZ3e39TUVPKLlfuiO8fKq6++GsuXL+/0MRobG2PUqFFl2R+AVRFEQEXtvPPOcfnll8eLL74Yzc3Nse2228Ypp5zS7hdpki4vv/xy7L///rH66qvHV77yldhjjz2qvUv0wFZbbdXuF5kWmzJlStxxxx1lea7uHCuf/OQn48477+z0MVb1O6AAysl3iACgzt19993tTlcttsYaa3R52fpyu//++2Px4sWd3j948OAuL3sOUE6CCAAAyCwXVQAAADJLEAEAAJkliAAAgMwSRAAAQGYJIgAAILP8HiKAOjZv3rxq70JdmDt3brV3AYB+4rLbACkmaOqDoAJIL0EEUCHihp4QUQCVIYgA+kDkkAbiCaD3BBFAJ8QO9UQ0AXRMEAGZJXiglWACskoQAXVJ7ED5iSagHgkioGaJHkgPsQTUKkEEpJrogdonloA0E0RAKggfyB6hBKSBIAIqSvgAqyKUgEoSREC/ET9AuYgkoL8IIqDPhA9QLUIJ6CtBBPSI+AHSTiQBPSGIgE6JH6BeiCSgM4IIiAjxA2SPSAIiBBFklgACKCWQIJsEEWSA+AHoHZEE9U8QQZ0RPwD9SyRBfRFEUOMEEEB1CSSobYIIaowAAkg3gQS1RRBBygkggNomkCDdBBGkjAACqG8CCdJFEEGVCSCAbBNIUF2CCCpMAAHQFYEElSWIoAJEEAC9IY6g/wki6AcCCID+IJCg/AQRlIkIAqCSxBGUhyCCPhBBAKSBOILeE0TQAwIIgFogkKD7BBGsgggCoJaJI+iaIIIOiCAA6pE4gvYEEbxHBAGQJeIIVmqo9g4AAABUixkiMs2sEACYLSLbBBGZI4IAoHPiiKwRRGSCCAKAnhNHZIEgom6JIAAoH3FEvRJE1A0BBACVI5CoF4KImieEAKB6hBG1ThBRk0QQAKSPOKIWCSJqhggCgNohjqgVgojUE0IAULuEEWkniEglEQQA9UcckUaCiFQRQgBQ/4QRaSKIqDoRBADZJY6oNkFE1QghACBPGFEtgoiKEkEAwKqIIypJEFERQggA6ClhRCUIIvqVEAIA+koY0Z8EEWUnggCA/iKOKDdBRNkIIQCgUoQR5SKI6DMhBABUizCirwQRvSaEAIC0EEb0liCix4QQAJBWwoieEkR0mxACAGqFMKK7BBFdEkEAQK0TR3RFENEhIQQA1BthREcEESWEEABQ74QRxQQRESGEAIDsEUZECKLME0IAQNYJo2wTRBklhAAASgmjbGqo9g5QeWIIAKA9Y6RsMkOUIX7IAQC6x2xRdgiiDBBCAAC9I4zqnyCqY0IIAKA8hFH9EkR1SAgBAPQPYVR/BFEdEUIAAJUhjOqHIKoDQggAoPJEUX0QRDVMCAEApIM4ql2CqAYJIQCAdBJGtUcQ1RAhBABQG4RR7Wio9g7QPWIIAKB2GLvVDjNEKeeHCQCgtpktSjdBlFJCCACgvgijdHLKXAqJIQCA+mOMl05miFLEDwkAQDaYLUoPM0QpIYYAALLD2C89zBBVmR8GAIBsM1tUXWaIqkgMAQBgTFhdZoiqwEEPAEBHzBZVnhmiChNDAAB0xlix8swQVYiDGwCAnjBbVBlmiCpADAEA0FPGkJVhhqgfOYgBACgHs0X9xwxRPxFDAACUi7Fl/zFDVGYOVgAA+pPZovIyQ1RGYggAgP5mzFlegqhMHJgAAFSKsWf5OGWujxyMAABUk1Po+sYMUR+IIQAAqs2YtG8EUS858AAASAtj095zylwPOdgAAEgzp9D1jBmiHhBDAACknTFrzwiibnJgAQBQK4xdu08QdYMDCgCAWmMM2z2+Q9QFBxEAAPXA94o6Z4aoE2IIAIB6YWzbOUHUAQcMAAD1xhi3Y4KoDQcKAAD1yli3PUEEAABkliAqopgBAKh3xrylBNF7HBgAAGSFsW+rzF9228EAAECWZf2S3JmeIRJDAABkXdbHxJkNoqx/8AAAkJflsXEmgyjLHzgAAHQkq2PkzAVRVj9oAABYlSyOlTMVRFn8gAEAoCeyNmbOTBBl7YMFAIDeytLYORNBlKUPFAAAyiErY+i6D6KsfJAAAFBuWRhL13UQZeEDBACA/lTvY+q6DaJ6/+AAAKBS6nlsXZdBVM8fGAAAVEO9jrHrLojq9YMCAIBqq8exdl0FUT1+QAAAkCb1NuaumyCqtw8GAADSqp7G3nURRPX0gQAAQC2olzF4zQdRvXwQAABQa+phLF7TQVQPHwAAANSyWh+T13QQAQAA9EUuSZKk2jvRU7VeoQAAUI/mzp1b7V3osZqbIRJDAACQTrU4Vq+pIKrFNxgAALKk1sbsNRNEtfbGAgBAVtXS2L0mgqiW3lAAAKB2xvA1EUQAAAD9IfVBVCtlCQAAlKqFsXyqg6gW3kAAAKBzaR/TpzaI0v7GAQAA3ZPmsX1qgwgAAKC/pTKI0lyQAABAz6V1jJ+6IErrGwUAAPRNGsf6qQqiNL5BAABA+aRtzJ+qIAIAAKik1ARR2koRAADoH2ka+6ciiNL0hgAAAP0vLQ2QiiACAACohqoHUVrKEAAAqKw0tEBVgygNbwAAAFA91W6Cqs8QAQAAVEvVgqjaJQgAAKRDNdugKkEkhgAAgGLVagSnzAEAAJlV8SAyOwQAAHSkGq1ghggAAMisigaR2SEAAKArlW4GM0QAAEBmVSyIzA4BAADdUcl2MEMEAABkVkWCyOwQAADQE5VqiH4PIjEEAAD0RiVawilzAABAZgkiAAAgs/o1iJwuBwAA9EV/N4UZIgAAILP6LYjMDgEAAOXQn21hhggAAMisfgkis0MAAEA59VdjmCECAAAyq+xBZHYIAADoD/3RGmaIAACAzBJEAABAZpU1iJwuBwAA9KdyN4cZIgAAILMEEQAAkFllCyKnywEAAJVQzvYwQwQAAGSWIAIAADJLEAEAAJlVliDy/SEAAKCSytUgZogAAIDMEkQAAEBm9TmInC4HAABUQzlaxAwRAACQWYIIAADILEEEAABkVp+CyPeHAACAauprk5ghAgAAMksQAQAAmSWIAACAzOp1EPn+EAAAkAZ9aRMzRAAAQGYJIgAAILMEEQAAkFmCCAAAyCxBBAAAZFavgsgV5gAAgDTpbaOYIQIAADJLEAEAAJkliAAAgMwSRAAAQGYJIgAAILMEEQAAkFmCCAAAyCxBBAAAZJYgAgAAMksQAQAAmSWIAACAzBJEAABAZgkiAAAgswQRAACQWYIIAADILEEEAABkliACAAAySxABAACZJYgAAIDMEkQAAEBmCSIAACCzBBEAAJBZgggAAMgsQQQAAGSWIAIAADJLEAEAAJkliAAAgMzqVRDNnTu33PsBAADQa71tFDNEAABAZgkiAAAgswQRAACQWYIIAADILEEEAABkVq+DyJXmAACANOhLm5ghAgAAMksQAQAAmSWIAACAzOpTEPkeEQAAUE19bRIzRAAAQGYJIgAAILMEEQAAkFl9DiLfIwIAAKqhHC1ihggAAMgsQQQAAGRWWYLIaXMAAEAllatBzBABAACZJYgAAIDMEkQAAEBmlS2IfI8IAACohHK2hxkiAAAgswQRAACQWWUNIqfNAQAA/anczWGGCAAAyCxBBAAAZFbZg8hpcwAAQH/oj9YwQwQAAGRWvwSRWSIAAKCc+qsxzBABAACZ1W9BZJYIAAAoh/5sCzNEAABAZvVrEJklAgAA+qK/m8IMEQAAkFmCCAAAyKx+DyKnzQEAAL1RiZaoyAyRKAIAAHqiUg3hlDkAACCzKhZEZokAAIDuqGQ7mCECAAAyq6JBZJYIAADoSqWbwQwRAACQWRUPIrNEAABAR6rRCmaIAACAzKpKEJklAgAAilWrEao2QySKAACAiOq2gVPmAACAzKpqEJklAgCAbKt2E1R9hqjabwAAAFAdaWiBqgcRAABAtaQiiNJQhgAAQOWkpQFSEUQR6XlDAACA/pWmsX9qgggAAKDSUhVEaSpFAACg/NI25k9VEEWk7w0CAADKI41j/dQFUUQ63ygAAKD30jrGT2UQAQAAVEJqgyitBQkAAPRMmsf2qQ2iiHS/cQAAwKqlfUyf6iCKSP8bCAAAdKwWxvKpDyIAAID+UhNBVAtlCQAAtKqVMXxNBFFE7byhAACQdbU0dq+ZIIqorTcWAACyqNbG7DUVRBG19wYDAEBW1OJYPZckSVLtneitefPmVXsXAAAg82oxhPJqboYIAACgXGo6iGq5RAEAoB7U+pi8poMoovY/AAAAqFX1MBav+SCKqI8PAgAAakm9jMHrIogi6ucDAQCAtKunsXfdBFFEfX0wAACQRvU25q6rIIqovw8IAADSoh7H2nUXRBH1+UEBAEA11esYuy6DKKJ+PzAAAKi0eh5b120QRdT3BwcAAJVQ72Pqug6iiPr/AAEAoL9kYSxd90EUkY0PEgAAyikrY+hMBFFEdj5QAADoqyyNnTMTRBHZ+mABAKA3sjZmzlQQRWTvAwYAgO7K4lg5c0EUkc0PGgAAupLVMXImgygiux84AAC0leWxcWaDKCLbHzwAAEQYE+eSJEmqvRNpMG/evGrvAgAAVEzWQygv0zNExRwQAABkhbFvK0EEAABkliAqopQBAKh3xrylBFEbDhAAAOqVsW57gqgDDhQAAOqNMW7HBFEnHDAAANQLY9vOuex2N7gkNwAAtUgIrZoZom5wIAEAUGuMYbtHEHWTAwoAgFph7Np9gqgHHFgAAKSdMWvP+A5RL/leEQAAaSKEescMUS854AAASAtj094TRH3gwAMAoNqMSfvGKXNl4hQ6AKrt0tWHRkTE7NferPKeAJUghMrDDFGZOCABqKZ8DAHZYOxZPoKojByYAFSDGIJsMeYsL6fM9ROn0AFQCZ3FkNPmoP4Iof5hhqifOGAB6G9mhiA7jC37jxmiCjBbBEC5rSqGzBBBfRBC/c8MUQU4kAEop+7MDJk9gtpnDFkZZogqzGwRAH3Rk9AxSwS1SQhVlhmiCnOAA9BbZn2g/hkrVp4ZoioyWwRAd/QlhMwSQW0QQtVjhqiKHPgArIpZIah/xoTVZYYoJcwWAdBWOWLIDBGklxBKBzNEKeEHAoBi5ZoZMsME6WTslx5miFLIbBFAtpU7YswSQXoIofQxQ5RCflAAsqs/ZnTMEkE6GOOlkxmilDNbBJAN/R0tZomgeoRQugmiGiGMAOpXJWZwBBFUnhCqDU6ZqxF+oADqU6VOZ3PaHFSWsVvtMENUg8wWAdSHSkeKWSLof0Ko9giiGiaMAGpTNWdrRBH0DyFUuwRRHRBGALWj2qeuCSIoHxFUHwRRHRFGAOlW7RjKE0XQd2KofgiiOiSMANIlLSGUJ4ig94RQ/RFEdUwYAVRf2mIoQhBBbwih+iWIMkAYAVRHGmMoTxRB9wih+ieIMkQYAVRGmkMoTxBB14RQdgiiDBJGAP2nFmIoTxRBe0IoexqqvQNUnh90gPK7dPWhNRVDQHvGSNlkhijjzBYB9F0th5BZIhBCWSeIiAhhBNAbtRxCeYKILBNCRAgi2hBGAN1TDzGUJ4rIGiFEMUFEh4QRQMfqKYTyBBFZIYToiCCiS8IIoFU9xlCeKKKeCSG6IojoNnEEZFU9h1CeIKLeiCC6SxDRY8IIyJIsxFCeKKIeCCF6ShDRa8IIqGdZCqFioohaJYToLUFEnwkjoJ5kNYTyBBG1RgjRV4KIshFGQK3LegzliSJqgRCiXAQRZSeMgFoigtoTRKSZEKLcBBH9ShwBaSWEVk0YkRYiiP4kiKgIYQSkhRDqGVFENQkhKkEQUVHCCKgmMdRzgohqEEJUkiCiasQRUClCqG9EEZUggqgWQUTVCSOgvwih8hFF9BchRLUJIlJFHAHlIITKTxBRTiKINBFEpJIwAnpDCPUvUURfCSHSSBCReuIIWBUhVDmiiJ4SQaSdIKJmCCOgLSFUeYKI7hJC1ApBRE0SR5BtQqi6RBGdEUHUIkFEzRNHkB1CKD1EEXkiiFoniKgbwgjqlxBKJ1GUbUKIeiGIqFsCCWqbCKoNoig7BBD1ShCRCeIIaocQqi2CqL6JILJAEJE54gjSSQjVLlFUX0QQWSOIyDRxBNUnhOqDKKptIogsE0TwHnEE1SGIapsQql0iCFZqqPYOAAAAVIsZIuiA2SKoHDNEtcvsUO0xKwTtCSJYBXEE/UsQ1R4hVFtEEHRNEEEPiCMoP0FUO4RQ7RBB0H2CCPpAIEHfCaLaIIbSTQBB7wkiKBNxBL0jiNJNCKWXCILyEETQD8QRdJ8gSh8RlF4iCMpPEEEFCCTonCBKFzGULgII+p8gggoTR1BKEKWDEEoPEQSVJYigygQSWSeIqksIVZ8AguoSRJAyAomsEUTVIYSqRwBBuggiSDmBRL0TRJUlhCpPAEG6CSKoMQKJeiOIKkMIVY4AgtoiiKDGCSRqnSDqX0Ko/wkgqG2CCOqMQKLWCKL+IYT6jwCC+iKIIANEEmkmiMpLCJWX+IH6J4ggo0QSaSGIykMI9Z34gWwSREBECCSqRxD1ngjqGwEERAgioAsiiUoQRD0nhHpO/ACdEURAj4gkyk0QdZ8Q6h7xA/SEIAL6TCTRF4KoayKoa+IH6CtBBPQboUR3CKL2RFB7wgfoL4IIqCiRRFuCaCUR1Er8AJUkiIBUEErZleUgynoECR8gDQQRkGpCqf5lLYiyGEHCB0gzQQTULLFUH7IQRFmIINED1CpBBNQlsVQ76jGI6jWARA9QjwQRkFmiKR3qJYhqPYLEDpBVggigE4KpMmoxiGoxfgQPQMcEEUAfiKa+q4UgSnsAiR2A3hNEABUinjqWtiBKS/yIHIDKEEQAKZaFiKp0EFUjeMQNQHoJIoA6VgtBVe4g6o/gETQA9UsQAQAAmdVQ7R0AAACoFkEEAABkliACAAAySxABAACZJYgAAIDMEkQAAEBmCSIAACCzBBEAAJBZgggAAMgsQQQAAGSWIAIAADJLEAEAAJkliAAAgMwSRAAAQGYJIgAAILMEEQAAkFmCCAAAyCxBBAAAZJYgAgAAMksQAQAAmSWIAACAzBJEAABAZgkiAAAgswQRAACQWYIIAADIrAHV3gGoR8uWLYvly5dXezcAqDNNTU0xaNCgau8G1BVBBGW2bNmyGDF4jVgey6q9KwDUmbFjx8azzz4riqCMBBGU2fLly2N5LIuPxa4xINccuYbcyjtyDUW33/t3Q65wO9fQULS8oXW9/P25hpXrF2+fy7VZN4qWt1237fad70uSy7WeUFv8XF3e/97y926XrNt2WUPx/UWP897t4sdPSpa3rptE/nYUXlfr/W3WLVpeWK/4+fPLGzrevqBk+05ud/BelDx/h+t2cDs6ub/NvnS6fWfLVvH8eZ0u6+CxOnovIpd0e18it/JVFV5bu3WTjp+zeHnh+VuX5braPpKiH4HWZ851sn1+ecljFm2fa7N90Y9LNJRs33p/Q9GyhvdeffHj5NdtaHM7IqIhSpc1dHA7/1id3Z9/zpXLWlqfK9re3xKNRdu0rrvy8RsjiVzR9q3rFi0rvv3euvnnacy1FB6z8b3nyz9u4bk6eKzGXEthHxuL1sv/GDdG8ePmt0kKj7Vy+/x2rY/TWPT6G4v2Jf+5FB4rWt/L/P0rl0Xre5Xfl1xE43ufSOuyXDQUlrXebszllzUULVt5+/U3WmLcFgti+fLlggjKSBBBPxkQA2NAbmDkioKk+PbKf7eOmnK5oiAqvr+h6P5VBlGu3XadBlGuzf19DqJcawSUBE9rGJQ1iIoHxv0cRB1v38ntwva5ottF268iSFYVEf0eRB3cn9c2iLp6X3oVREXPVfEgKr7dwfZtg6jtNq2HS+vAuVdB1ME2vQ2i0uDpfhAVL1/5786CqDgCeh5EDR0GUdLJ7e4HUWNhv3LR8N4bmo+RlUGUv50ripSkaFlS9FhR2JfW529d1lUQNfYgiBq7FUS+9g39xU8XAACQWYIIAADILEEEAABkliACAAAySxABAACZJYgAAIDMEkQAAEBmCSIAACCzBBEAAJBZgggAAMgsQQQAAGSWIAIAADJrQLV3AOrVu/FORNIQuST33pLi2+/9O8kVbueShqLlDa3rtbx3f67o/lzxvxtabxcePtfBum23b3N/kivcTnK5iCTabL+q+99bHhHREqXr5l9yfllD8f1Fj1PYlVzhZSUly1vXzT995CKioe1jtVm3w7e9aJsO3vaS/YqOtu/kdtv3Ile6vON1O7gdndzfZl863b6zZat4/rxOl3XwWB29F5FLur0vkVv5qgqvrd26ScfPWby88Pyty3JdbR9J0Y9A6zPnOtk+v7zkMYu2z7XZvujHJZKS7VvvT4qWJe+9+vzjtBRt09DmdkREQ5Qua+jgduG/PJ3c3xDFy1panyva3t8SjUXbtK678vEbI4lc0fat6xYtK7793rr552nMtRQes/G958s/buG5OnisxlxLYR8bi9bL/xg3RvHj5rdJCo+1cvv8dq2P01j0+huL9iX/WRUeK1rfy/z9K5dF63uV35dcRON7n0jrslw0FJa13m7Mta7Xumzl473+RksA5SeIoMySJImhQ4fGXW/euHKEt6LaewRAvRg6dGgkSbLqFYFuE0RQZrlcLt588834xz/+EcOHD6/27gBQJ15//fVYd911I1c8cw30mSCCfjJ8+HBBBACQci6qAAAAZJYgAgAAMksQQZk1NzfH3Llzo7m5udq7AkAd8f8X6B+5xKVKAACAjDJDBAAAZJYgAgAAMksQAQAAmSWIAACAzBJEAABAZgkiKKPvf//7MX78+Bg0aFBsvfXW8X//93/V3iUAatypp54aW221VQwbNixGjx4de+65Zzz++OPV3i2oG4IIyuSKK66IL33pSzF37tx44IEHYtNNN41p06bFSy+9VO1dA6CG3XnnnTFnzpy477774ne/+1288847sdNOO8Vbb71V7V2DuuD3EEGZbL311rHVVlvFeeedFxERLS0tse6668aRRx4ZJ554YpX3DoB68fLLL8fo0aPjzjvvjE984hPV3h2oeWaIoAyWL18e999/f0ydOrWwrKGhIaZOnRr33ntvFfcMgHqzZMmSiIhYc801q7wnUB8EEZTBokWLYsWKFTFmzJiS5WPGjIkXX3yxSnsFQL1paWmJY445JiZPnhyTJk2q9u5AXRhQ7R0AAKB75syZE4888kjcdddd1d4VqBuCCMpg5MiR0djYGAsXLixZvnDhwhg7dmyV9gqAenLEEUfEDTfcEL///e/j/e9/f7V3B+qGU+agDJqammKLLbaIW2+9tbCspaUlbr311th2222ruGcA1LokSeKII46Ia6+9Nm677bZYb731qr1LUFfMEEGZfOlLX4pZs2bFlltuGR/96Efje9/7Xrz11ltx0EEHVXvXAKhhc+bMiZ///Ofxq1/9KoYNG1b4buqIESNi8ODBVd47qH0uuw1ldN5558UZZ5wRL774Ymy22WZxzjnnxNZbb13t3QKghuVyuQ6XX3LJJTF79uzK7gzUIUEEAABklu8QAQAAmSWIAACAzBJEAABAZgkiAAAgswQRAACQWYIIAADILEEEAABkliACAAAySxABAACZJYgAAIDMEkQAAEBm/X/xDb1RF5UeTAAAAABJRU5ErkJggg==" + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "## Filter catalog and plot filtered pixels\n", + "\n", + "filtered_catalog = catalog.filter_by_cone(ra, dec, radius)\n", + "\n", + "inspection.plot_pixels(filtered_catalog)" + ], + "metadata": { + "collapsed": false + } } ], "metadata": { diff --git a/src/.pylintrc b/src/.pylintrc index 387e4914..5ff594bc 100644 --- a/src/.pylintrc +++ b/src/.pylintrc @@ -183,6 +183,7 @@ good-names=i, k, ex, Run, + ra, _ # Good variable names regexes, separated by a comma. If names match any regex, diff --git a/src/hipscat/catalog/association_catalog/association_catalog.py b/src/hipscat/catalog/association_catalog/association_catalog.py index 4c0a4e54..29ef7f77 100644 --- a/src/hipscat/catalog/association_catalog/association_catalog.py +++ b/src/hipscat/catalog/association_catalog/association_catalog.py @@ -2,11 +2,9 @@ import pandas as pd -from hipscat.catalog.catalog_type import CatalogType -from hipscat.catalog.association_catalog.association_catalog_info import ( - AssociationCatalogInfo, -) +from hipscat.catalog.association_catalog.association_catalog_info import AssociationCatalogInfo from hipscat.catalog.association_catalog.partition_join_info import PartitionJoinInfo +from hipscat.catalog.catalog_type import CatalogType from hipscat.catalog.dataset.dataset import Dataset from hipscat.io import FilePointer, paths @@ -59,8 +57,6 @@ def _read_args( cls, catalog_base_dir: FilePointer ) -> Tuple[CatalogInfoClass, JoinPixelInputTypes]: args = super()._read_args(catalog_base_dir) - partition_join_info_file = paths.get_partition_join_info_pointer( - catalog_base_dir - ) + partition_join_info_file = paths.get_partition_join_info_pointer(catalog_base_dir) partition_join_info = PartitionJoinInfo.read_from_file(partition_join_info_file) return args + (partition_join_info,) diff --git a/src/hipscat/catalog/catalog.py b/src/hipscat/catalog/catalog.py index 197d3d9e..06074da1 100644 --- a/src/hipscat/catalog/catalog.py +++ b/src/hipscat/catalog/catalog.py @@ -1,6 +1,7 @@ """Container class to hold catalog metadata and partition iteration""" from __future__ import annotations +import dataclasses from typing import Tuple, Union import pandas as pd @@ -10,6 +11,7 @@ from hipscat.catalog.dataset.dataset import Dataset from hipscat.catalog.partition_info import PartitionInfo from hipscat.io import FilePointer, file_io, paths +from hipscat.pixel_math.cone_filter import filter_pixels_by_cone from hipscat.pixel_tree.pixel_tree import PixelTree from hipscat.pixel_tree.pixel_tree_builder import PixelTreeBuilder @@ -82,9 +84,7 @@ def get_pixels(self): return self.partition_info.data_frame @classmethod - def _read_args( - cls, catalog_base_dir: FilePointer - ) -> Tuple[CatalogInfoClass, PartitionInfo]: + def _read_args(cls, catalog_base_dir: FilePointer) -> Tuple[CatalogInfoClass, PartitionInfo]: args = super()._read_args(catalog_base_dir) partition_info_file = paths.get_partition_info_pointer(catalog_base_dir) partition_info = PartitionInfo.read_from_file(partition_info_file) @@ -96,3 +96,21 @@ def _check_files_exist(cls, catalog_base_dir: FilePointer): partition_info_file = paths.get_partition_info_pointer(catalog_base_dir) if not file_io.does_file_or_directory_exist(partition_info_file): raise FileNotFoundError(f"No partition info found where expected: {str(partition_info_file)}") + + def filter_by_cone(self, ra: float, dec: float, radius: float) -> Catalog: + """Filter the pixels in the catalog to only include the pixels that overlap with a cone + + Args: + ra (float): Right Ascension of the center of the cone in degrees + dec (float): Declination of the center of the cone in degrees + radius (float): Radius of the cone in degrees + + Returns: + A new catalog with only the pixels that overlap with the specified cone + """ + filtered_cone_pixels = filter_pixels_by_cone(self.pixel_tree, ra, dec, radius) + filtered_catalog_info = dataclasses.replace( + self.catalog_info, + total_rows=None, + ) + return Catalog(filtered_catalog_info, filtered_cone_pixels) diff --git a/src/hipscat/pixel_math/cone_filter.py b/src/hipscat/pixel_math/cone_filter.py new file mode 100644 index 00000000..996ab2d0 --- /dev/null +++ b/src/hipscat/pixel_math/cone_filter.py @@ -0,0 +1,51 @@ +import healpy as hp +import numpy as np +import pandas as pd + +from hipscat.catalog.partition_info import PartitionInfo +from hipscat.pixel_tree import PixelAlignment, PixelAlignmentType, align_trees +from hipscat.pixel_tree.pixel_tree import PixelTree +from hipscat.pixel_tree.pixel_tree_builder import PixelTreeBuilder + + +def filter_pixels_by_cone(pixel_tree: PixelTree, ra: float, dec: float, radius: float) -> PixelTree: + """Filter the leaf pixels in a pixel tree to return a partition_info dataframe with the pixels + that overlap with a cone + + Args: + ra (float): Right Ascension of the center of the cone in degrees + dec (float): Declination of the center of the cone in degrees + radius (float): Radius of the cone in degrees + + Returns: + A catalog_info dataframe with only the pixels that overlap with the specified cone + """ + max_order = max(pixel_tree.pixels.keys()) + cone_tree = _generate_cone_pixel_tree(ra, dec, radius, max_order) + cone_alignment = align_trees(pixel_tree, cone_tree, alignment_type=PixelAlignmentType.INNER) + pixels_df = cone_alignment.pixel_mapping[ + [PixelAlignment.PRIMARY_ORDER_COLUMN_NAME, PixelAlignment.PRIMARY_PIXEL_COLUMN_NAME] + ] + filtered_pixels_df = pixels_df.drop_duplicates() + partition_info_df = filtered_pixels_df.rename( + columns={ + PixelAlignment.PRIMARY_ORDER_COLUMN_NAME: PartitionInfo.METADATA_ORDER_COLUMN_NAME, + PixelAlignment.PRIMARY_PIXEL_COLUMN_NAME: PartitionInfo.METADATA_PIXEL_COLUMN_NAME, + } + ) + return partition_info_df.reset_index(drop=True) + + +def _generate_cone_pixel_tree(ra: float, dec: float, radius: float, order: int): + """Generates a pixel_tree filled with leaf nodes at a given order that overlap with a cone""" + n_side = hp.order2nside(order) + center_vec = hp.ang2vec(ra, dec, lonlat=True) + radius_radians = np.radians(radius) + cone_pixels = hp.query_disc(n_side, center_vec, radius_radians, inclusive=True, nest=True) + cone_pixel_info_dict = { + PartitionInfo.METADATA_ORDER_COLUMN_NAME: np.full(len(cone_pixels), order), + PartitionInfo.METADATA_PIXEL_COLUMN_NAME: cone_pixels, + } + cone_partition_info_df = pd.DataFrame.from_dict(cone_pixel_info_dict) + cone_tree = PixelTreeBuilder.from_partition_info_df(cone_partition_info_df) + return cone_tree diff --git a/src/hipscat/pixel_tree/pixel_alignment.py b/src/hipscat/pixel_tree/pixel_alignment.py index 0d91afb1..add1a058 100644 --- a/src/hipscat/pixel_tree/pixel_alignment.py +++ b/src/hipscat/pixel_tree/pixel_alignment.py @@ -8,7 +8,6 @@ from hipscat.pixel_tree.pixel_tree import PixelTree from hipscat.pixel_tree.pixel_tree_builder import PixelTreeBuilder - LEFT_TREE_KEY = "left" RIGHT_TREE_KEY = "right" @@ -84,9 +83,7 @@ def align_trees( """ tree_builder = PixelTreeBuilder() - pixels_to_search = _get_children_pixels_from_trees( - [left, right], left.root_pixel.pixel - ) + pixels_to_search = _get_children_pixels_from_trees([left, right], left.root_pixel.pixel) while len(pixels_to_search) > 0: search_pixel = pixels_to_search.pop(0) @@ -96,14 +93,10 @@ def align_trees( if left_node.node_type == right_node.node_type: if left_node.node_type == PixelNodeType.LEAF: # Matching leaf nodes get added to the aligned tree - tree_builder.create_node_and_parent_if_not_exist( - search_pixel, PixelNodeType.LEAF - ) + tree_builder.create_node_and_parent_if_not_exist(search_pixel, PixelNodeType.LEAF) else: # For matching inner nodes search into their children to check for alignment - pixels_to_search += _get_children_pixels_from_trees( - [left, right], search_pixel - ) + pixels_to_search += _get_children_pixels_from_trees([left, right], search_pixel) else: # Nodes with non-matching types: one must be a leaf and the other an inner node if left_node.node_type == PixelNodeType.LEAF: @@ -112,45 +105,33 @@ def align_trees( else: tree_with_leaf_node = RIGHT_TREE_KEY inner_node = left_node - if _should_include_all_pixels_from_tree( - tree_with_leaf_node, alignment_type - ): + if _should_include_all_pixels_from_tree(tree_with_leaf_node, alignment_type): # If the alignment type means fully covering the tree with the leaf node, then # create a leaf node in the aligned tree and split it to match the partitioning # of the other tree to ensure the node is fully covered - tree_builder.create_node_and_parent_if_not_exist( - search_pixel, PixelNodeType.LEAF - ) + tree_builder.create_node_and_parent_if_not_exist(search_pixel, PixelNodeType.LEAF) tree_builder.split_leaf_to_match_partitioning(inner_node) else: # Otherwise just add the subtree from the inner node to include all the # overlapping pixels - tree_builder.create_node_and_parent_if_not_exist( - search_pixel, PixelNodeType.INNER - ) + tree_builder.create_node_and_parent_if_not_exist(search_pixel, PixelNodeType.INNER) tree_builder.add_all_descendants_from_node(inner_node) elif search_pixel in left and search_pixel not in right: # For nodes that only exist in one tree, include them if the alignment type means that # tree should have all its nodes included if _should_include_all_pixels_from_tree(LEFT_TREE_KEY, alignment_type): - tree_builder.create_node_and_parent_if_not_exist( - search_pixel, left[search_pixel].node_type - ) + tree_builder.create_node_and_parent_if_not_exist(search_pixel, left[search_pixel].node_type) tree_builder.add_all_descendants_from_node(left[search_pixel]) elif search_pixel in right and search_pixel not in left: if _should_include_all_pixels_from_tree(RIGHT_TREE_KEY, alignment_type): - tree_builder.create_node_and_parent_if_not_exist( - search_pixel, right[search_pixel].node_type - ) + tree_builder.create_node_and_parent_if_not_exist(search_pixel, right[search_pixel].node_type) tree_builder.add_all_descendants_from_node(right[search_pixel]) tree = tree_builder.build() pixel_mapping = _generate_pixel_mapping_from_tree(left, right, tree) return PixelAlignment(tree, pixel_mapping, alignment_type) -def _get_children_pixels_from_trees( - trees: List[PixelTree], pixel: HealpixInputTypes -) -> List[HealpixPixel]: +def _get_children_pixels_from_trees(trees: List[PixelTree], pixel: HealpixInputTypes) -> List[HealpixPixel]: """Returns the combined HEALPix pixels that have child nodes of the given pixel from trees This returns a list of HEALPix pixels, not the actual child nodes, and does not contain @@ -173,9 +154,7 @@ def _get_children_pixels_from_trees( return list(pixels_to_add) -def _should_include_all_pixels_from_tree( - tree_type: str, alignment_type: PixelAlignmentType -) -> bool: +def _should_include_all_pixels_from_tree(tree_type: str, alignment_type: PixelAlignmentType) -> bool: """If for a given alignment type, the left or right tree should include all pixels or just the ones that overlap with the other tree. @@ -188,13 +167,12 @@ def _should_include_all_pixels_from_tree( """ left_add_types = [PixelAlignmentType.OUTER, PixelAlignmentType.LEFT] right_add_types = [PixelAlignmentType.OUTER, PixelAlignmentType.RIGHT] - return (tree_type == LEFT_TREE_KEY and alignment_type in left_add_types) or \ - (tree_type == RIGHT_TREE_KEY and alignment_type in right_add_types) + return (tree_type == LEFT_TREE_KEY and alignment_type in left_add_types) or ( + tree_type == RIGHT_TREE_KEY and alignment_type in right_add_types + ) -def _generate_pixel_mapping_from_tree( - left: PixelTree, right: PixelTree, aligned: PixelTree -) -> pd.DataFrame: +def _generate_pixel_mapping_from_tree(left: PixelTree, right: PixelTree, aligned: PixelTree) -> pd.DataFrame: """Generates a pixel mapping dataframe from two trees and their aligned tree The pixel mapping dataframe contains columns for the order and pixel of overlapping pixels in @@ -226,27 +204,15 @@ def _generate_pixel_mapping_from_tree( right_leaf_nodes = [None] for left_node in left_leaf_nodes: for right_node in right_leaf_nodes: - pixel_mapping_dict[PixelAlignment.ALIGNED_ORDER_COLUMN_NAME].append( - leaf_node.hp_order - ) - pixel_mapping_dict[PixelAlignment.ALIGNED_PIXEL_COLUMN_NAME].append( - leaf_node.hp_pixel - ) + pixel_mapping_dict[PixelAlignment.ALIGNED_ORDER_COLUMN_NAME].append(leaf_node.hp_order) + pixel_mapping_dict[PixelAlignment.ALIGNED_PIXEL_COLUMN_NAME].append(leaf_node.hp_pixel) left_order = left_node.hp_order if left_node is not None else None left_pixel = left_node.hp_pixel if left_node is not None else None - pixel_mapping_dict[PixelAlignment.PRIMARY_ORDER_COLUMN_NAME].append( - left_order - ) - pixel_mapping_dict[PixelAlignment.PRIMARY_PIXEL_COLUMN_NAME].append( - left_pixel - ) + pixel_mapping_dict[PixelAlignment.PRIMARY_ORDER_COLUMN_NAME].append(left_order) + pixel_mapping_dict[PixelAlignment.PRIMARY_PIXEL_COLUMN_NAME].append(left_pixel) right_order = right_node.hp_order if right_node is not None else None right_pixel = right_node.hp_pixel if right_node is not None else None - pixel_mapping_dict[PixelAlignment.JOIN_ORDER_COLUMN_NAME].append( - right_order - ) - pixel_mapping_dict[PixelAlignment.JOIN_PIXEL_COLUMN_NAME].append( - right_pixel - ) + pixel_mapping_dict[PixelAlignment.JOIN_ORDER_COLUMN_NAME].append(right_order) + pixel_mapping_dict[PixelAlignment.JOIN_PIXEL_COLUMN_NAME].append(right_pixel) pixel_mapping = pd.DataFrame.from_dict(pixel_mapping_dict) return pixel_mapping diff --git a/src/hipscat/pixel_tree/pixel_tree.py b/src/hipscat/pixel_tree/pixel_tree.py index 9afa674f..5ad53901 100644 --- a/src/hipscat/pixel_tree/pixel_tree.py +++ b/src/hipscat/pixel_tree/pixel_tree.py @@ -78,9 +78,7 @@ def get_node(self, pixel: HealpixInputTypes) -> PixelNode | None: def __getitem__(self, item): return self.get_node(item) - def get_leaf_nodes_at_healpix_pixel( - self, pixel: HealpixInputTypes - ) -> List[PixelNode]: + def get_leaf_nodes_at_healpix_pixel(self, pixel: HealpixInputTypes) -> List[PixelNode]: """Lookup all leaf nodes that contain or are within a given HEALPix pixel - Exact matches will return a list with only the matching pixel @@ -111,9 +109,7 @@ def get_leaf_nodes_at_healpix_pixel( return [] return [node_in_tree] - def _find_first_lower_order_leaf_node_in_tree( - self, pixel: HealpixInputTypes - ) -> PixelNode | None: + def _find_first_lower_order_leaf_node_in_tree(self, pixel: HealpixInputTypes) -> PixelNode | None: pixel = get_healpix_pixel(pixel) for delta_order in range(1, pixel.order + 1): lower_pixel = pixel.convert_to_lower_order(delta_order) diff --git a/src/hipscat/pixel_tree/pixel_tree_builder.py b/src/hipscat/pixel_tree/pixel_tree_builder.py index 0cb68607..e76c6add 100644 --- a/src/hipscat/pixel_tree/pixel_tree_builder.py +++ b/src/hipscat/pixel_tree/pixel_tree_builder.py @@ -6,10 +6,7 @@ from hipscat.catalog.partition_info import PartitionInfo from hipscat.pixel_math import HealpixPixel -from hipscat.pixel_math.healpix_pixel_convertor import ( - HealpixInputTypes, - get_healpix_pixel, -) +from hipscat.pixel_math.healpix_pixel_convertor import HealpixInputTypes, get_healpix_pixel from hipscat.pixel_tree.pixel_node import PixelNode from hipscat.pixel_tree.pixel_node_type import PixelNodeType from hipscat.pixel_tree.pixel_tree import PixelTree @@ -174,9 +171,7 @@ def create_node( node_to_replace = None if pixel in self: if not replace_existing_node: - raise ValueError( - f"Cannot create node at {str(pixel)}, node already exists" - ) + raise ValueError(f"Cannot create node at {str(pixel)}, node already exists") node_to_replace = self[pixel] parent.remove_child_link(node_to_replace) node = PixelNode(pixel, node_type, parent) @@ -192,7 +187,6 @@ def create_node( for child in node_to_replace.children: self._remove_node_and_children_from_tree(child.pixel) - def remove_node(self, pixel: HealpixInputTypes): """Remove node in tree @@ -263,9 +257,5 @@ def split_leaf_to_match_partitioning(self, node_to_match: PixelNode): parent_node = self[node.parent.pixel] if parent_node.node_type == PixelNodeType.LEAF: parent_node.node_type = PixelNodeType.INNER - for child_pixel in parent_node.pixel.convert_to_higher_order( - delta_order=1 - ): - self.create_node( - child_pixel, PixelNodeType.LEAF, parent=parent_node - ) + for child_pixel in parent_node.pixel.convert_to_higher_order(delta_order=1): + self.create_node(child_pixel, PixelNodeType.LEAF, parent=parent_node) diff --git a/tests/hipscat/catalog/test_catalog.py b/tests/hipscat/catalog/test_catalog.py index 6794cdaf..3fa842f6 100644 --- a/tests/hipscat/catalog/test_catalog.py +++ b/tests/hipscat/catalog/test_catalog.py @@ -2,6 +2,7 @@ import os +import pandas as pd import pandas.testing import pytest @@ -73,6 +74,45 @@ def test_load_catalog_small_sky_order1(small_sky_order1_dir): assert len(cat.get_pixels()) == 4 +def test_cone_filter(small_sky_order1_catalog): + filtered_catalog = small_sky_order1_catalog.filter_by_cone(315, -66.443, 0.1) + assert len(filtered_catalog.partition_info.data_frame) == 1 + assert filtered_catalog.partition_info.data_frame[PartitionInfo.METADATA_PIXEL_COLUMN_NAME][0] == 44 + assert filtered_catalog.partition_info.data_frame[PartitionInfo.METADATA_ORDER_COLUMN_NAME][0] == 1 + assert (1, 44) in filtered_catalog.pixel_tree + assert len(filtered_catalog.pixel_tree.pixels[1]) == 1 + assert filtered_catalog.catalog_info.total_rows is None + + +def test_cone_filter_big(small_sky_order1_catalog): + filtered_catalog = small_sky_order1_catalog.filter_by_cone(315, -66.443, 30) + assert len(filtered_catalog.partition_info.data_frame) == 4 + assert (1, 44) in filtered_catalog.pixel_tree + assert (1, 45) in filtered_catalog.pixel_tree + assert (1, 46) in filtered_catalog.pixel_tree + assert (1, 47) in filtered_catalog.pixel_tree + + +def test_cone_filter_multiple_order(catalog_info): + partition_info_df = pd.DataFrame.from_dict( + { + PartitionInfo.METADATA_ORDER_COLUMN_NAME: [6, 7, 7], + PartitionInfo.METADATA_PIXEL_COLUMN_NAME: [30, 124, 5000], + } + ) + catalog = Catalog(catalog_info, partition_info_df) + filtered_catalog = catalog.filter_by_cone(47.1, 6, 30) + assert len(filtered_catalog.partition_info.data_frame) == 2 + assert (6, 30) in filtered_catalog.pixel_tree + assert (7, 124) in filtered_catalog.pixel_tree + + +def test_cone_filter_empty(small_sky_order1_catalog): + filtered_catalog = small_sky_order1_catalog.filter_by_cone(0, 0, 0.1) + assert len(filtered_catalog.partition_info.data_frame) == 0 + assert len(filtered_catalog.pixel_tree) == 1 + + def test_empty_directory(tmp_path): """Test loading empty or incomplete data""" ## Path doesn't exist