{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# pyphot - A tool for computing photometry from spectra\n", "\n", "\n", "Some examples are provided in this notebook\n", "\n", "Full documentation available at http://mfouesneau.github.io/docs/pyphot/" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "%matplotlib inline\n", "\n", "import pylab as plt\n", "import numpy as np\n", "\n", "import sys\n", "sys.path.append('../')\n", "from pyphot import sandbox as pyphot\n", "from pyphot import Vega, Sun" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## The Sun and Vega" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "`pyphot` provides convenient interfaces to a spectral representation of the Sun and Vega.\n", "In this notebook, we show how they can be used." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Library contains: 256 filters\n" ] } ], "source": [ "# get the internal default library of passbands filters\n", "lib = pyphot.get_library()\n", "print(\"Library contains: \", len(lib), \" filters\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Vega" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Suppose one has a calibrated spectrum and wants to compute the vega magnitude throug the HST WFC3 F110W passband," ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Vega magnitude of Vega in HST_WFC3_F110W is : 0.000000 mag\n", "AB magnitude of Vega in HST_WFC3_F110W is : 0.751950 mag\n", "ST magnitude of Vega in HST_WFC3_F110W is : 2.372749 mag\n" ] } ], "source": [ "# convert to magnitudes\n", "import numpy as np\n", "\n", "# We'll use Vega spectrum as example\n", "vega = Vega()\n", "f = lib['HST_WFC3_F110W']\n", "# compute the integrated flux through the filter f\n", "# note that it work on many spectra at once\n", "fluxes = f.get_flux(vega.wavelength, vega.flux, axis=-1)\n", "# convert to vega magnitudes\n", "mags = -2.5 * np.log10(fluxes.magnitude) - f.Vega_zero_mag\n", "print(\"Vega magnitude of Vega in {0:s} is : {1:f} mag\".format(f.name, mags))\n", "mags = -2.5 * np.log10(fluxes.magnitude) - f.AB_zero_mag\n", "print(\"AB magnitude of Vega in {0:s} is : {1:f} mag\".format(f.name, mags))\n", "mags = -2.5 * np.log10(fluxes.magnitude) - f.ST_zero_mag\n", "print(\"ST magnitude of Vega in {0:s} is : {1:f} mag\".format(f.name, mags))" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(,\n", " ,\n", " )" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "f.AB_zero_Jy, f.Vega_zero_Jy, f.ST_zero_Jy" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## The Sun" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The internal reference to the Solar spectrum comes in two flavors: an observed one and a theoretical one.\n", "By default, the interface is set to theoretical.\n", "\n", "In addition, the Sun is at $1\\,au$ but can be set to any distance. Below we instanciate the two flavors and also the Sun if it was at $10\\, pc$ (absolute flux units)." ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": [ "from pyphot import unit\n", "sun_obs = Sun(flavor='observed')\n", "sun_th = Sun() # default is theoric spectrum\n", "sun_th_10pc = Sun(distance=10 * unit['pc'])" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "WARNING: VerifyWarning: Invalid keyword for column 2: Column null option (TNULLn) must be an integer for binary table columns (got 1.6e+38). The invalid value will be ignored for the purpose of formatting the data in this column. [astropy.io.fits.column]\n", "WARNING: VerifyWarning: Invalid keyword for column 3: Column null option (TNULLn) must be an integer for binary table columns (got 1.6e+38). The invalid value will be ignored for the purpose of formatting the data in this column. [astropy.io.fits.column]\n", "WARNING: VerifyWarning: Invalid keyword for column 4: Column null option (TNULLn) must be an integer for binary table columns (got 1.6e+38). The invalid value will be ignored for the purpose of formatting the data in this column. [astropy.io.fits.column]\n", "WARNING: VerifyWarning: Invalid keyword for column 5: Column null option (TNULLn) must be an integer for binary table columns (got 1.6e+38). The invalid value will be ignored for the purpose of formatting the data in this column. [astropy.io.fits.column]\n", "WARNING: VerifyWarning: Invalid keyword for column 7: Column null option (TNULLn) must be an integer for binary table columns (got 1.6e+38). The invalid value will be ignored for the purpose of formatting the data in this column. [astropy.io.fits.column]\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY4AAAEOCAYAAACetPCkAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xd4VNXWwOHfSiPUAKEpRSJShYgYQDoKAiooICKKigKi2K/XiheFz16uBQEBpVkuCCoKyBUUxNCUckVaQEOPtNACAdL398eZhCGmzEwyOTOT9T7PPDOz58w5K0fMyj57n7XFGINSSinlqiC7A1BKKeVfNHEopZRyiyYOpZRSbtHEoZRSyi2aOJRSSrlFE4dSSim3aOJQSinlFk0cSiml3KKJQymllFs0cSillHJLiN0BeEO1atVM/fr17Q5DKaX8yoYNG44aY6oXtl1AJo769euzfv16u8NQSim/IiJ7XdlOL1UppZRyiyYOpZRSbtHEoZRSyi0BOcahlPJf6enpJCQkkJKSYncoASs8PJw6deoQGhrq0fc1cSilfEpCQgIVK1akfv36iIjd4QQcYwzHjh0jISGBqKgoj/ahl6qUUj4lJSWFyMhITRpeIiJERkYWqUeniUNd4M8/t7N7dzyZ6WmcOH7M7nBUKaVJw7uKen71UlUplpVlOHUujQNJqURVK8/AFyeyoMy/cj6v4rTtr1lNaBu0ne8y25BepRHlz+zjbOP+dOs9iArlypZ88EqVoD179tC7d2+2bNlidygX6Nq1K2+//TYxMTElelxNHKXYhI8n88iBZ6jseL+gzPnPfslqytVBcTnv2wZtB+DG4LVwaq3VuC0Wtj0OQFxWPTZFDaNTnyFcXM055Sil8pKRkUFIiH/+CtZLVaXQseRUTp5N45EDz/zts68yO8GYJFqOXsWjaQ/xQugTTL3oRQCSTTivVn2Zz6o9xpLMq/hf1mWszWoMQNOgfdy290UuHl+f02MuZst7/dix4ksyMzNL9GdTqji88847NG/enObNm/Pee+8B1i/6IUOGEB0dzYABAzh79iwAzz77LM2aNSM6Oponn3wSgMTERG655RZat25N69atWbVqFQBjxoxhxIgR9OjRg7vvvpu2bduydevWnON27dqVDRs2cObMGYYOHUrr1q258sor+fbbbwE4d+4cgwYNIjo6mttuu41z586V5GnJ4Z/pTnlk98GjjPpgGlU5zbXB/+OWYIjLqsufPT7lph+7AnDj6G8ACA8N5t8vvUyQCMFBAjxBBWCU0/6SzqVToUwIR5KSCD+9n0rTOgKQFlSW5ieXwdJlsHQYO4IbkhLRgOo3v8zFlzQs0Z9Z+bexC7ay7cCpYt1ns4sr8WKfy/P9fMOGDUyfPp1ff/0VYwxt27alS5cu7Nixg6lTp9KhQweGDh3KxIkTGTp0KPPmzWP79u2ICCdPngTgscce4x//+AcdO3Zk37599OzZk7i4uJz9r1y5krJly/Luu+8yZ84cxo4dy8GDBzlw4ABXXXUVo0aN4tprr2XatGmcPHmSNm3a0L17dyZPnky5cuXYtGkTmzZtolWrVsV6blyliSOALf9lHefKXUz1MhnEzIomCpgVduE2a8PbM6TjlfCj9T487Pw/idDggjukEWWtOeA1qlSGKpVJGfI9qccPEHnVLezZ/SfHYj+mxl8/UD91N2WO/wnTv2dLVn1WlO9Bh74jaN7wMoKCdBBU+ZaVK1fSr18/ypcvD0D//v1ZsWIFdevWpUOHDgDceeedjBs3jscff5zw8HCGDx/OjTfeSO/evQH48ccf2bZtW84+T506xenTpwG46aabKFvWGhccOHAg1113HWPHjmXOnDnceuutACxZsoT58+fz9ttvA9ZMs3379hEbG8ujjz4KQHR0NNHR0SVwRv5OE0eAOn5oL12/7w7A3IzOxOTzX7oh+wH4qdt8IitXoSj/DMOj2hHumBZeP6oh9aPeAN4gKzOLLctnk7Z7NVX2LWHkuSkwawoAGYRwZvBCIhq2K8KRVaAqqGfgLcaYPNtzz0QSEUJCQli7di1Lly5l9uzZjB8/nmXLlpGVlcWaNWtyEoSz7IQEULt2bSIjI9m0aRNffPEFkydPzonhq6++onHjxoXGYQcd4whAmfE/8cn4MTnvbw2JzXfbSo4ZUdd06kJ0C+/89RIUHETzbnfQavh46o/Zzo7e83I+CyGDiM97wZgIUl6qzaG9O7wSg1Ku6ty5M9988w1nz57lzJkzzJs3j06dOrFv3z7WrFkDwKxZs+jYsSPJyckkJSVxww038N5777Fx40YAevTowfjx43P2md2el0GDBvHmm2+SlJREixYtAOjZsycffPBBThL77bffcmL7/PPPAdiyZQubNm0q/hPgAk0cgeK3z2HnMjK+e4rgz/ryeMjXF3w8NO1JVmc2Y23Hjy9ob9Dp1pKMEgkKonHMtWwctpdZnZZwsMVIzomVvMIzk6k1vQ2/vtCW5V+8x/GTSSUam1IArVq14p577qFNmza0bduW4cOHU6VKFZo2bcrMmTOJjo7m+PHjjBw5ktOnT9O7d2+io6Pp0qUL7777LgDjxo1j/fr1REdH06xZMyZNmpTv8QYMGMDs2bMZOHBgTtvo0aNJT08nOjqa5s2bM3r0aABGjhxJcnIy0dHRvPnmm7Rp08a7JyMfkl+3zJ/FxMSYUrcex5iIPJuzWgwkqHoj9l3+EGMXbGXc7VeS9V40Fc/9BUDmCycdg9/2MhlpxP84lR3xf9DsyCIuDTrESVOeuZldiOgwlIG9uoMPdNGV98XFxdG0aVO7wwh4eZ1nEdlgjCn0phAd4wgETsk/wwTRM+0NOldL5umOVSjb9l4A6gFT72ltbfTYL0z97yo+X5vAMh9IGgASEkbDXiNpCCSnvMna1QtJ/3Uq96Z8T8ivi4hbU4/YsE606/8Q0c1K/rq3Uuo87XH4qwWPwcZZMPoImQd+J3hKZ/7Mqs3AtNGcoBKLH+9M41oV7Y6yyF77fBFnti1hUPBPNA/aA8DqzGbQehgxPQYTFq53rQca7XGUDO1xlEYbZljPR+MxC58AIMFU4wSVWPd8d6pXLJP/d/3Ic4NvAG7AGMPLH31Gx/0fcVnQX9T53z/hf/9kekZP6vZ4hG6dOvrEbBOlSgNNHP7uw3Yk1e1BJJDU5nF+aNM5YJKGMxHhXyPuAu7CZGYwe+rbDDrwGveGLIZli9nwY0N217uFXreNpELFyoXuTynlOZ1V5e8y04jcsxCAKhc3omFN/788VRgJDuHW4c+w+q6d/HHXBl5Jv4PKJDMg4XUq/PsS1r/QmuWxP+U7H18pVTTa4wggnZtebHcIJSY4SGjfoBpQjedf+ZCszEy+/vglmv71JTFBf8Cyvmz/sS5TMm7k7Zdf1zvUlSpG2uPwUwnBdf/WJsFheWxZOgQFB9P//jE0fHETj0a8z6SM3jQJ2s87YZOIfbELr0+ZaXeIyk+cPHmSiRMnArB8+fKcMiIlZfny5axevTrn/aRJk/jkk0882leFChWKK6wLaOLwU6Em9e+NpThxZAsJDmLcP+7hgZc/Z3rTj/gyszOtg7bz7IFHmT+6F+u2xdsdovJxzonDWzIyMvL9LHfieOCBB7j77ru9Go+7NHH4qTJZeSz7GOzZwvOB6t7bBjLgpQUk3r+J7zLbcFPwGlp8cTULJjxJ0tk8Eq9SWGXSd+7cScuWLXnqqadITk5mwIABNGnShMGDB+eMnW3YsIEuXbpw1VVX0bNnTw4ePAhY5UWuvvpqoqOj6devHydOnACskumjRo2iS5cuvP/++3mWXt+zZw+TJk3i3XffpWXLlqxYsYIxY8bkFDuMj4+ne/fuXHHFFbRq1YqdO3eSnJxMt27daNWqFS1atMgpwe5NPj/GISKXAs8DEcaYAXbHY7vTh+Hfjchz3pBOR81T/YtrUf+lH3h36kxa7JlOn8SPOPzGXFY2f4Ebbx1md3iqIP99Fg5tLt591moB17+e78evv/46W7ZsYePGjSxfvpybb76ZrVu3cvHFF9OhQwdWrVpF27ZteeSRR/j222+pXr06X3zxBc8//zzTpk3j7rvv5oMPPqBLly688MILjB07NmdNj5MnT/Lzzz8DcMcdd+RZev2BBx6gQoUKOWt7LF26NCe2wYMH8+yzz9KvXz9SUlLIysoiLCyMefPmUalSJY4ePcrVV1/NTTfd5NXp6bYkDhGZBvQGjhhjmju19wLeB4KBj40xrxtjdgHDRORLO2L1OX9tyHmZRihhpNsYjH/5x7AhpGfcydMvPsOboR9x49YnSNjyKnv7zadDS70bXeWtTZs21KlTB4CWLVuyZ88eKleuzJYtW7juuusAyMzM5KKLLiIpKYmTJ0/SpUsXAIYMGZJTKh3gtttuy3ldUOn1vJw+fZq//vqLfv36ARAeHg5Aeno6o0aNIjY2lqCgIP766y8OHz5MrVq1iukM/J1dPY4ZwHggZ8RHRIKBCcB1QAKwTkTmG2O25bmH0mj2YDh8frWwlDKRhKUesjEg/xMaEsybr7zNgOcu4qnQObQN2k6VedcwackQ7nzkJSqUDbc7ROWsgJ5BSSlT5vx9UcHBwWRkZGCM4fLLL8+plpstKangwpzOJdULKr2el/yml3/++eckJiayYcMGQkNDqV+/PikpeVzKLka2jHEYY2KB47ma2wDxxphdxpg0YDZwc4kH58u2L4QTu3PextUZZGMw/u29p0dS5cEfGJb2TzZkNeKBs1PY/1pr3vlYZ1+VdhUrVizwL3+Axo0bk5iYmJM40tPT2bp1KxEREVSpUoUVK1YA8Omnn+b0PnLLr/R6fsevVKkSderU4ZtvrFU6U1NTOXv2LElJSdSoUYPQ0FB++ukn9u7d6/4P7SZfGhyvDY5VhSwJQG0RiRSRScCVIvJcfl8WkREisl5E1icmJno71pJx5ihs+xZ+HAvHd/3t46rhWSSFVLMhMP9Xp0o5GtWqxNRXX6DDmFieDn6KSnKGJxIeZeXo9sTt9P7/fMo3RUZG0qFDB5o3b85TTz2V5zZhYWF8+eWXPPPMM1xxxRW0bNkyZybUzJkzeeqpp4iOjmbjxo288MILee4jv9Lrffr0Yd68eTmD484+/fRTxo0bR3R0NO3bt+fQoUMMHjyY9evXExMTw+eff06TJk2K8WzkzbYihyJSH1iYPcYhIrcCPY0xwx3v7wLaGGMecXffAVPkcPoNsNda5D6l9UOEr5twwcfm+rdgxTtI8kFoeSf0nZDXXpSLjp84wZZ3+tA5eDPHTEXGhQ5jzKgXkSBf+vsq8GmRw5JRlCKHvvR/RALgfFdbHeCATbHYLzkxJ2kAf0sah658HGk97PxEqo6Pl2BwgalqlSqYu+bxSr2P2WtqMjbjPZa9eC2/bd5id2hK+RRfShzrgIYiEiUiYcAgYL7NMdkjPQXevqzATaTdwxAUDKmOa6HlIksgsMDXpVF1nh96K3X/uYI35R7aBW3jsi+7M+r5J0gv4KYtpUoTWxKHiMwC1gCNRSRBRIYZYzKAh4HFQBwwxxiztaD9BKx1HxW6SWSkI1F0c1w/DdeKsMWpekQ5nn7xfd6ImsbvWZfyauhUtr/Snl3b1tkdmlK204WcfE3aGXjVhWKFY3Q97pJy6lwaUz54laFnP6Yqp/ml8o20GD6Z8hUCvxKxHeLi4mjSpImur+JFxhi2b98eEGMcCuCHFy94O5seNgWislUqG8aTT4/BjFzNruAorj75HSnj2sLOZXaHFpDCw8M5duyYlsX3EmMMx44dy7mB0BM+X3Kk1Ek/B8D9aY9z4pJezLm/HYyJsDkoBRBZsx6Rozfyx5oFXLr2Rfi0H0QPgp6vQHmdFl1c6tSpQ0JCAgEzrd4HhYeH59wN7wlNHD4mlWBOm0qcufQGPhlSaI9R2aBRuz4Qcx2seBtWvgfxP0CvN6DFAK0XVgxCQ0OJioqyOwxVAL1U5WNSUtPIIJhbY+oQHhpsNT7uKPLWsKd9gakLhYbDtf+C+2OhShR8PRxmDYKkv+yOTCmv08ThY1LT0skkiOoVnNYNr1wP/vkH3PYp1GgGkQ3tC1BdqGYzGLYEer4Gu36GiVfD+umg1+dVANPE4WNS0zMwCNUqlrnwg4o1IaQMPLgGHvHTGWOBKigY2j1o/be5uCUsfBxm9smzTIxSgUATh49JSUsnywjVKpQpfGPlW6pGwd3zoc/7cPB3mNge1kyArEy7I1OqWGni8DHJ59IgKJgq5XQ1P78kAlfdAw/+Apd2gcWjYFpPOLLd7siUKjaaOHxIRmYWJ86kEBISojc/+buI2nD7bOj/MRzbCZM7wa+TIVPLlij/p4nDh/z8RyJnU9OpEB5mdyiqOIhA9K3w0FqI6gz/fRo+usa6jKWUH9PE4UO+23SQsCBDhbKaOAJKheow+Eu4dQacPgRTroHYtyAj1e7IlPKIJg4fkZyawbot2+kUspWgCM/v6FQ+SgQu7wcPr4XL+8Kyl2FSR9i/1u7IlHKbJg4PJadm8MfhgpeXdMfz8zbTJvM3ymYmn694qwJP2Spwy1S4Y45VXmZqD2vsQyk/oonDQ8NnrqPHu7HFUojt8KkUvt+UwL/DrKUjqXF5kfepfJgINOpp3ffR8Dpr7GPB4zl1ypTydZo4PPTLruMApGZkFXlfK/48ysigr883BGsJsVKhTEUY9B9oNQQ2TIcPO8DeNXZHpVShNHF4KHu2bGp60RJHZpbhybm/83iII3E0H1DEyJRfCQ6Fm8bBnV9BVjpM7wXzH4Uzx+yOTKl8aeLwUJAjc6RmFO2u4KVxhwnBmttvyleHAVOLHJvyQ5d1t24abPcw/PYZTGgNm+ZqzSvlkzRxeCg4J3EUrcexZNthuobHAyDXv1HkuJQfCytvre3xwMrzFXf/cxskJdgdmVIX0MThoSDHmUtJ97zHkZVl+PmPRIZEbLQaarYohsiU33OuuLtnBUy4GtZ+BFlFH09Tqjho4vBQUXscZ9MyiB67hIjknXRKmm81lqlQXOEpf+dccbdua1j0pDX+cXib3ZEppYnDU0FBVuLwtMexKSGJ5NQMugRtOt8YpolD5VKlPtz5NfSdBMfiYXJnWPYKpKfYHZkqxTRxeCioiD2OzQlJANzf+Oz5Rk0cKi8i0PJ2eGgdNL8FYt/Uu86VrXw+cYjIpSIyVUS+tDsWZ8FF6HEYY5i7YT91KodT/dDP0LAHDPvx/MCJUnkpHwn9J1tTd9PPWeXal/xLbxxUJc6rv6lEZJqIHBGRLbnae4nIDhGJF5FnC9qHMWaXMWaYN+P0RFD2fRwe9DgOJqXwx+FknrgiAzl7DJr2sa5jK+WKy7pbYx+t7obVH1iXr/avszsqVYp4+0/cGUAv5wYRCQYmANcDzYDbRaSZiLQQkYW5HjW8HJ/HshzT6zOy3J9n/9u+kwD0/3Wg1VC9SXGFpUqL8ErWSoN3zYO0szCtBywZrWMfqkR4NXEYY2KB47ma2wDxjp5EGjAbuNkYs9kY0zvX44g34yuKLMeNWVkeJY4TXBpy9HxD+erFFZYqbRpca/U+rrwLVo+zFozS3ofyMjsuqtcG9ju9T3C05UlEIkVkEnCliDxXwHYjRGS9iKxPTEwsvmjzkZ0wPOpx7D/JS+Wdhmw0caiiCK/kKFvy9fnexw8vaO9DeY0diSOvNVHz/e1rjDlmjHnAGNPAGPNaAdtNMcbEGGNiqlf3/i/i7IAz3bwpK+7gKX7ff5LaZc5AuUh4aqfev6GKx2Xdzvc+Vr1v9T4S1tsdlQpAdiSOBKCu0/s6wAEb4igS4+EYx1uLdxBRNpR6WQegYU8oX80L0alSK3fvY+p12vtQxc6OxLEOaCgiUSISBgwC5tsQR5Fkj3FkupE4jDGs3X2cV2ouIyj5IOhKf8pbLusGD66GK+/U3ocqdt6ejjsLWAM0FpEEERlmjMkAHgYWA3HAHGPMVm/G4Q3ZiSMj0/XEsXb3cZJTM4g5vdxqiL7NC5Ep5RAeATd9oL0PVey8umKQMeb2fNoXAYu8eWxvy+5ouNPjuG3KLwCElK0Ala+Gapd5IzSlLpTd+1jyL6v3seN76DsR6sTYHZnyU3qrsofcnVWVkXl+EL3SuQSoGuWVuJTKU07v4ytIS3b0Pl6EjFS7I1N+SBOHh86Pcbg2q2rPsTMAvNKnIUGnD1rrLShV0rLvOm85GFa9B5O7wIHf7I5K+RlNHB4wxrh95/h3mw4hAl2qnwUMVL3UewEqVZDwCLh5PAz+ElJOwkfdrIq7GWl2R6b8hCYODziPaxR25/ifh08zY9VuYv9MJLp2BHVSrdX+9FKVsl3D66zeR/RAq+LuR9fCoc12R6X8gCYODzj3MrJfn05J51jy368X3zn1V8Ys2MaGvSfo0rgGrHOsKa49DuULylaBfpNg0H8g+TBMuQZi34bMDLsjUz5ME4cH0pwGurN7H9e8vZyrXv7xb9umO03X7d6kGiRuh/qdoFxV7weqlKua3AgP/mJVal72klW2JPEPu6NSPkoThwfSnUqpZ/c4jibnfX24XFgwAPd3vpQWZQ7DueN6/4byTeUj4dbpMGA6HN9t3TS4ejxkebbKpQpcmjg84NyLyMwyfLxiV57bpWVkcTQ5leEdo3juhqbIYcd9jrVblUSYSnmmeX+r99HgWljyPMy4EY7n/W9clU6aODyQnunc48ji5e/i8tzuf/tOkJKeRZsox2WpzHTrObSct0NUqmgq1rTGPfp+CIe3wYcdYO1H4GZRTxWYNHF4IK8xjtyMMfxx+DQALepEOBodXf6gYK/Gp1SxEIGWd1gzr+pdDYuehM/6wcn9hX9XBTSvlhwJVBf0OHLVqnrlu200qF6BY2fSeGvxDmpXLkutSuHWh1mOmSpBetqVH4mobdW72jAdFv8LPmwPvV6zbiKUvFZJUIFOf4N5oFxoCD2a1WTJtsN/63F8tGI3AE1qVQTg/26+HMn+nyt7kFG0x6H8jAjEDIVLr4FvH7IecQus5Wsr1rI7OlXC9FKVB+pFlmPK3TFEVStPpsn7UlXCiXPc3e4SujWteb4xO3Foj0P5q6pRMGQh9HwNdi2HCW1h85fnF6hRpYImjiIoExLE2bS8pyomp2bQ/OKICxv/XGI9B+lpV34sKAjaPQgPrIRqDeGrYTB3CJw5andkqoQU+KeviLgybzTdGFMq6xTUqBTOkVP5r23Qs7lTFz4zA+J/sF5rj0MFgmoN4d7vYfU4WP4a7FllXbpq2tvuyJSXFfYb7GesFfsKGgGLAuoXV0D+5KJK4cT+kZjnZ+Nuv5KIsqHnG/bEnn8t2uNQASI4BDo9AY16wbz74YvB1g2u179hlTNRAamwxLHOGHNtQRuIyLJijMev1K5SNud1zUplOHzKqlX11oBobrri4gs3Tjtz/nVIWZQKKDWbwX3LrDpXsW/B7li4aTw07G53ZMoLCvzTt7Ck4eo2garZRZVyXt/T/ny12w6XVfv7xs533uoYhwpEwaFwzXNw31IIrwyf3wLzH4GUU3ZHpoqZyxfbRSQa65JUzneMMV97ISa/0a1pjZzX5csE89XIdnz9v7/O37fh7E/H+MaDv5RQdErZ5OIrYcRyWP4qrP4Adi6HvhMgqrPNgani4tKfviIyDZgG3AL0cTxK/QiYiNClUXUAyoWFcNUlVXmlXwuCgnINCZ3cD3tWwDXPQ42mNkSqVAkLDYfr/g+GLrZ6IjP7wKKnIe2s3ZGpYuBqj+NqY0wzr0bip6pVKAPAubQC1i/Ysch6bjGgBCJSyofUbWNN2106Fn6dBPE/WvWv6rW1OzJVBK5ebF8jIpo48nBjtDXl9qKIAga8962BiLq6eJMqncLKWbOshiywCn1O7wVLRkN6/lPZlW9zNXHMxEoeO0Rkk4hsFpFN3gwsm4g0FZFJIvKliIwsiWO649omNVnx9DUXjHf8TcIGqBNTckEp5YuiOsODq+HKu6x7P6Z0gQO/2R2V8oCriWMacBfQi/PjG30K+5KITBORIyKyJVd7L0cSiheRZwvahzEmzhjzADAQ8MnfvnWrljtfjyovyYetHodSpV2ZinDTOBj8FaQkwcfdYfkb55ccUH7B1cSxzxgz3xiz2xizN/vhwvdmYCWbHCISDEwArgeaAbeLSDMRaSEiC3M9aji+cxOwEljq6g/mM9LPQWYqlK1sdyRK+Y6G3a1y7Zf3t2ZfTb0OEnfYHZVykauJY7uI/EdEbheR/tmPwr5kjIkFjudqbgPEG2N2GWPSgNnAzcaYzcaY3rkeRxz7mW+MaQ8MduNn8w1Jf1nP4REFb6dUaVO2CtzyEdw6E07shcmdYc1EXSzKD7g6q6oskAr0cGozgCf3cdQGnFeCSQDynWIhIl2B/kAZYFEB240ARgDUq1fPg7C8ZP006zm0vL1xKOWrLu8L9drBgkdh8XPWLMSbJ0CVS+yOTOXDpcRhjLm3GI+Z12BAvjWZjTHLgeWF7dQYMwWYAhATE+M7NZ7THaVGLu9rbxxK+bKKNeH22fDbZ/D9c9ZStb1egyvv1MWifJBLiUNEwoFhwOVAzm3RxpihHhwzAXAeKa4DHPBgP/5hwwzrOVTrUylVIBFodZc1++qbB2H+w7B9IfQZZyUW5TNcHeP4FKgF9MSqmFsHOO3hMdcBDUUkSkTCgEHAfA/35ducCxsqpVxT5RLrno/sxaImXg1bv7E7KuXE1cRxmTFmNHDGGDMTuBFoUdiXRGQWsAZoLCIJIjLMGJMBPAwsBuKAOcaYrZ6F7+OyF7bpUuCMY6VUbtmLRd0fayWSuUPgq+Fw7oTdkSlcHxzPnmR9UkSaA4dwYQ0OY8zt+bQvooCB7oBx1pE4LrrC3jiU8lfVG8OwH2DFOxD7JuxZCTePh8u0XLudXO1xTBGRKsBorMtK24A3vRZVoEhNtp7DKxW8nVIqf8Gh0PUZGP6jNa39s1tg4T/O//+lSpxLicMY87Ex5oQx5mdjzKXGmBrGmEneDs7vZThq8ejCTUoV3cVXwoifod3DsH46TOoI+3SZAjsUtub4EwV9box5p3jDCTDpjhLSOqNKqeIRGg49X4HGN8A3D8C0XtDhUWvJgpAydkdXahTW46hYPDNoAAAWSUlEQVRYyEMVJP2c9ayJQ6niVb8DjFwNre6GVe/DlK5wsETqrioKHxwvZ4x5RkRuNcbMLZGIAkn2NdjQcvbGoVQgyi6Y2KS3dc/HR9daYyEd/gHBLi9uqjxQWI/jBhEJBZ4riWACzn+fsp61x6GU9zTqYS3J3LQPLHsZpvWEo/F2RxXQCksc3wNHgWgROeX0OC0iugJ9QZwLtZXRq3pKeVW5qnDrdLhlKhyLtwbOf52iBRO9pMDEYYx5yhgTAXxnjKnk9KhojNE5pgU5sdt6rt9Ja+0oVVJaDLB6H/U7Wj3+T/tCUoLdUQUcV6fj3uztQALO3lXWc89X7I1DqdKm0kUweC70eR8S1sPEdrDxP2B8p/apvyswcYjIwsJ24Mo2pdIfi6HyJVAr2u5IlCp9ROCqe2DkKqjZHL4ZCV/cCcmJdkcWEAqbetBRRAoqQChYq/gpZ8bA/rXQ4Fq9TKWUnapGwT0LYc0EWPaSVTCxz3vWQLryWGGJw5VLVGnFEUhASVgPZ45A3dZ2R6KUCgq2bhJseB3Mu9/qeVxxO/R6XZd09lCBicMY83NJBRJQvr7Peq6jiUMpn1GjKQxfCrFvQezbsDvWWmmwwTV2R+Z3XC1yqNwR2cB61qq4SvmW4FC4ZpRVcTe0nDXratFTkHbW7sj8iiaO4mYMHNoCzW+xOxKlVH7qXAUPrIC2I2HtFOu+j/3r7I7Kb7iUOESkRh5tjYs/nABwLB6SD1n3byilfFdoWbj+dWu1wcw0mNYDlv4fZOiwbWFc7XGsEJGB2W9E5J/APO+E5Od2LbeeL+1iaxhKKRdFdbYKJl5xB6z4t1Xz6tAWu6Pyaa4mjq7AXSIyV0RigUZAG69F5c92/wwRdaFKlN2RKKVcFV4J+k6AQbOsKwZTusLKdyEr0+7IfJKrd44fxKpb1Q5rydhPjDG6/FZu6Smw62e4tKvev6GUP2pyg1WypPH18OMYmH49HNtpd1Q+x9Uxjh+AtkBz4AbgXRF525uB+aWdSyH1FFze1+5IlFKeKl8NBn4C/T+CxO3WwPm6j7VkiRNXL1VNMMbcbYw5aYzZArQHkrwYl39a9zFUqAVROr6hlF8TgeiBMHIN1G0L3/0TPusPSX/ZHZlPcPVS1Te53mcYY17yTkh+Kukv2LkMYoZac8WVUv4vojbcNQ9u/Le1vvmH7WDTnFLf+3D1UtVpp7U4UkQkU0RKpMchIl1FZIWITBKRriVxTI/E/2A9aw0cpQKLCLQeDg+shOpNrMoQc4fAmWN2R2YbV3scFZ3W4ggHbgEmFPY9EZkmIkdEZEuu9l4iskNE4kXk2cIODyQD4YDvFtbf8T1E1LPKGiilAk9kA7j3v9B9DOz4r1Uwccd/7Y7KFh7dOe64dHWtC5vOAHo5N4hIMFbSuR6rsu7tItJMRFqIyMJcjxrACmPM9cAzwFhP4vW6E3usgfHG1+tsKqUCWVAwdPwH3PcTVKgBswbBNw9BSulaENWlFd1FpL/T2yAgBqsnUCBjTKyI1M/V3AaIN8bscux7NnCzMeY1oHcBuzsBlHEl3hIXv9S687TNfXZHopQqCbWaW8nj59et+z12x0LfiRBVOipGuNrj6OP06AmcxrWS63mpDex3ep/gaMuTiPQXkcnAp8D4ArYbISLrRWR9YmIJL9ayeS6UrwGRl5XscZVS9gkJg24vwNDF1oSYmb3h++cg/ZzdkXmdSz0OY8y9xXjMvK7l5Nt7McZ8DXxd2E6NMVOAKQAxMTElN+Xh7HHYtwa6PqeXqZQqjeq2sQom/jgGfpkI8T9Cv0lQ+yq7I/OaAhOHiHxAwb/UH/XgmAlAXaf3dYADHuzHN+xcZj3rvRtKlV5h5eGGt6DxDfDtQ/DxddDpn9Dl6YCcnl9Yj2O9F465DmgoIlHAX8Ag4A4vHKdkxC+FctWsm4SUUqVbg2usgonfPwuxb8Kfi6HfFKjRxO7IilVhYxzdjTEzgcrGmJm5H4XtXERmAWuAxiKSICLDjDEZwMPAYiAOmGOM2VrUH8QWmRnWoFjdNhCkS5sopbCWo+03CW77DJISYHJna83zrCy7Iys2hfU4rhKRS4ChIvIJucYnjDHHC/qyMeb2fNoXAYvcCdQn7V0JpxKg5yt2R6KU8jVN+1hXIuY/CotHWfd89J0IlevZHVmRFfZn8iSsqrhNgA25Ht64jOVf9q4GCYLLutkdiVLKF1WoAbfPstY2P7ARJraH3z73+5IlBSYOY8w4Y0xTYJox5lJjTJTT49ISitF3JW631t0oU9HuSJRSvkoErrwTRq6Ci6Lh2wfhizshuYRvGyhGrpYcGentQPzSsZ1QrZHdUSil/EGVS2DIQujxMvy5xCpZsv07u6PyiI7oFkXyYahYy+4olFL+IigI2j8CI36GShfB7Dv8smSJJg5PZWVZN/+Vi7Q7EqWUv6nZDIYvg05Pwu//gQ87wO4VdkflMk0cnkpNApOpiUMp5ZmQMOg2GoYugeAQR8mSUdYS1D5OE4enzjpmIperam8cSin/Vre1tdZH6+HwywSY0sWageXDNHF4KvW09Vymkr1xKKX8X1h5a5XBO7+ClCT4uBv8/JZ1k7EP0sThicPbrL8KAELL2huLUipwXNYdHlwDzfrCTy/DtJ5wNN7uqP5GE4cnpjmtTaWJQylVnMpWgQFTYcA0OBYPkzrC2o986qZBTRzuysywBsazaeJQSnlD81vgwV+gfgdY9CR82g+S/rI7KkATh/uO77zwfYgmDqWUl1S6CAZ/CTe+A/t/hQ/bwaa5tvc+NHG469DmC99rj0Mp5U0i0HqYNfOqWmP4ejh8ee/5mZ020MThrkObIMhpYZbQcvbFopQqPSIbwL3/tZarjVtolSz5Y4ktoWjicNfOn6BOzPn3oeH2xaKUKl2CQ6yVBe9bZt18/J9bYcFjkJpcomFo4nBHegoc3gL1O50vbqhjHEqpknZRNNz3E7R/FDbMhEkdYN8vJXZ4TRzuOL4LTBZUbwz3fAd3zdOV/5RS9ggNhx4vwb2LrMHy6dfDj2MgI9Xrh9bfeu44usN6rtbIWqClwbX2xqOUUpe0t9b6uPJOWPmutcCclxW2dKxydvRP6znyMnvjUEopZ2Uqwk0fQJsRUKuF1w+nPQ53HP0TIupCmM6kUkr5oBJIGqCJwz0n90HlS+yOQimlbKWJwx1JCVC5rt1RKKWUrTRxuCozHU4fgIg6dkeilFK28vnBcRHpBAzGirWZMaa9LYGcPmhNxY3QHodSqnTzao9DRKaJyBER2ZKrvZeI7BCReBF5tqB9GGNWGGMeABYCM70Zb4FO7reetcehlCrlvN3jmAGMBz7JbhCRYGACcB2QAKwTkflAMPBaru8PNcYccby+Axju5Xjzd8pRzlh7HEqpUs6ricMYEysi9XM1twHijTG7AERkNnCzMeY1oHde+xGRekCSMeZUfscSkRHACIB69eoVPfjcUh2HLlu5+PetlFJ+xI7B8drAfqf3CY62ggwDphe0gTFmijEmxhgTU7169SKGmIe0s9azllFXSpVydgyOSx5tBa5KYox50UuxuC79nPWsZdSVUqWcHT2OBMB5oKAOcMCGONyTfhaCy0BQsN2RKKWUrexIHOuAhiISJSJhwCBgvg1xuCf9nF6mUkopvD8ddxawBmgsIgkiMswYkwE8DCwG4oA5xpit3oyjWKSf0ctUSimF92dV3Z5P+yJgkTePXezSz2lxQ6WUQkuOuE4vVSmlFKCJw3VpeqlKKaVAE4frTu6DirXsjkIppWynicMVKafgxO4SWyRFKaV8mSYOVxx2TPqqFW1vHEop5QM0cbjisKO4b83m9sahlFI+QBOHKw5tgrJVodLFdkeilFK208ThikNboFZzkLzKbCmlVOmiiaMwmRlwZJuObyillIMmjsIci4eMFJ1RpZRSDpo4CqMD40opdQFNHIU5tAmCw6BaI7sjUUopn6CJozCHtkD1xhASZnckSinlEzRxFObQZh0YV0opJ5o4CnL6MJw5ogPjSinlRBNHQQ5vtp51YFwppXJo4ijIIUfiqKWJQymlsmniKMihLRBRF8pWsTsSpZTyGZo4CnJos45vKKVULpo48pN+Do79qYlDKaVy0cSRnyPbwGTpwLhSSuWiiSM/R+Ks55qX2xuHUkr5GJ9PHCLSTETmiMiHIjKgxA58JA5CwqFK/RI7pFJK+QOvJg4RmSYiR0RkS672XiKyQ0TiReTZQnZzPfCBMWYkcLfXgs3tyDar1EhQcIkdUiml/EGIl/c/AxgPfJLdICLBwATgOiABWCci84Fg4LVc3x8KfAq8KCI3AZFejve8I9shqnOJHU4ppfyFVxOHMSZWROrnam4DxBtjdgGIyGzgZmPMa0DvfHb1kCPhfO2tWC9w7iScPgA1mpbI4ZRSyp94u8eRl9rAfqf3CUDb/DZ2JJ5RQHngrQK2GwGMAKhXr17RIkzcbj1r4lBKqb+xI3HktXC3yW9jY8weHAmhIMaYKcAUgJiYmHz355Ij26xnTRxKKfU3dsyqSgDqOr2vAxywIY78HYmDsApWuRGllFIXsCNxrAMaikiUiIQBg4D5NsSRvyNxUL0JSF6dI6WUKt28PR13FrAGaCwiCSIyzBiTATwMLAbigDnGmK3ejMNtR+L0MpVSSuXD27Oqbs+nfRGwyJvH9lhyIpw9qolDKaXy4fN3jpe4REepEU0cSimVJ00cuR1xTMWt3sTeOJRSykdp4sjt6A4oUwkqXmR3JEop5ZM0ceSWuAOqNdIZVUoplQ9NHLkd/UMvUymlVAE0cTg7dwKSD0P1RnZHopRSPksTh7PEP6znao3tjUMppXyYJg5nR3dYz9U1cSilVH40cThL3GGt+le5iNV1lVIqgGnicGYM1I7RVf+UUqoAdpRV9129XrU7AqWU8nna41BKKeUWTRxKKaXcoolDKaWUWzRxKKWUcosmDqWUUm7RxKGUUsotmjiUUkq5RROHUkopt4gxxu4Yip2IJAJ7nZoigKRcm7nSVg04WuwBFiyvuLz5fVe2L2gbdz/T8+769nrei//7et7zjiHbJcaY6oV+2xgT8A9giidtwHpfiNWb33dl+4K2cfczPe963vW8+/95Ly2XqhYUoa2kFTUGd7/vyvYFbePuZ3reXd9ez3vxf1/Pu6VIMQTkpariIiLrjTExdsdR2uh5t4eed3v443kvLT0OT02xO4BSSs+7PfS828Pvzrv2OJRSSrlFexxKKaXcoolDKaWUWzRxKKWUcosmDjeISFMRmSQiX4rISLvjKU1EpLyIbBCR3nbHUlqISFcRWeH4N9/V7nhKCxEJEpFXROQDERlidzx5KfWJQ0SmicgREdmSq72XiOwQkXgReRbAGBNnjHkAGAj41fQ5X+POeXd4BphTslEGHjfPuwGSgXAgoaRjDSRunvebgdpAOj563kt94gBmAL2cG0QkGJgAXA80A24XkWaOz24CVgJLSzbMgDMDF8+7iHQHtgGHSzrIADQD1/+9rzDGXI+VtMeWcJyBZgaun/fGwBpjzBOAT17ZKPWJwxgTCxzP1dwGiDfG7DLGpAGzsf4KwBgz3xjTHhhcspEGFjfP+zXA1cAdwH0iUur/3XrKnfNujMlyfH4CKFOCYQYcN/+9J2Cdc4DMkovSdSF2B+CjagP7nd4nAG0d13n7Y/1PtMiGuAJdnufdGPMwgIjcAxx1+oWmikd+/977Az2BysB4OwILcHmed+B94AMR6QTE2hFYYTRx5E3yaDPGmOXA8pINpVTJ87znvDBmRsmFUqrk9+/9a+Drkg6mFMnvvJ8FhpV0MO7QLn/eEoC6Tu/rAAdsiqU00fNuDz3v9vDb866JI2/rgIYiEiUiYcAgYL7NMZUGet7toefdHn573kt94hCRWcAaoLGIJIjIMGNMBvAwsBiIA+YYY7baGWeg0fNuDz3v9gi0865FDpVSSrml1Pc4lFJKuUcTh1JKKbdo4lBKKeUWTRxKKaXcoolDKaWUWzRxKKWUcosmDlXqiMhyESnWsvgiUllEHnR631VEFroYyw5H1WWvE5G+2ZWevbT/n0QkubjPr/ItmjiUKh6VgQcL3Spvg40xJXXHcF+sEt5/IyJFrl1njLkGWF/U/SjfpolD+QQReVpEHnW8fldEljledxORzxyvPxSR9SKyVUTGOtquF5E5TvvpKiILHK97iMgaEfmfiMwVkQp5HDfPbURkj4iMdbRvFpEmjvbqIvKDo32yiOwVkWrA60ADEdkoIm85dl9BrNUit4vI5yKSV1G73PHcJyLrROR3EflKRMo52meIyDgRWS0iu0RkgKM9SEQmOs7JQhFZ5PTZ6yKyTUQ2icjbItIeuAl4yxFnA0eP51UR+Rl4TEQuEZGlju8sFZF6Tsf/0NGj2CUiXcRanChORGa4+Z9b+TtjjD70YfsDa72NuY7XK4C1QCjwInC/o72q4zkYq0pxNFaF531AecdnHwJ3AtWwSlJntz8DvOB4vRxrBceCttkDPOJ4/SDwseP1eOA5x+teWNV7qwH1gS1OP09XIAmrcF0QVrmJjnn83MuBGKf3kU6vX3aKYQYw17GvZljrOAAMwCrxHwTUwlrHYQBQFdjB+eoQlZ32MyDX8Sc6vV8ADHG8Hgp84/S92VgVXW8GTgEtHMfdALTM72fSR+A9tMehfMUG4CoRqQikYv2ijQE6YSUSgIEi8j/gN+ByoJmx6v18D/RxXGq5EfgWKxE1A1aJyEZgCHBJrmMWtk12SfENWIkBoCPWL1CMMd9zfsGdvKw1xiQYa/2QjU77KEhzsdb53oy1WNjlTp99Y4zJMsZsA2o6xTPX0X4I+MnRfgpIAT4Wa12NswUc8wun1+2A/zhef+rYf7YFxhgDbAYOG2M2O362rS7+bCpA6HocyicYY9JFZA9wL7Aa2IS18l8DIE5EooAngdbGmBOOyyPhjq9/ATyEtcLaOmPMacdloR+MMbcXcNjCtkl1PGdy/v+VQi835fH93PsoyAygrzHmd7EWruqaz/4k1/MFjDEZItIG6IZVdfVh4Np8jnmmgHici9llHz8rVyxZ6O+SUkV7HMqXxGIlh1isXsYDwEbHX7mVsH7BJYlITax1mrMtB1oB93H+r+dfgA4ichmAiJQTkUa5jufKNrmtBAY6tu8BVHG0nwYquvXT5q0icFBEQnFteeKVwC2OsY6aOBKNY6wmwhizCHgcaOlinKuxEg2O4690+ydQAU8Th/IlK4CLgDXGmMNYl1pWABhjfse6RLUVmAasyv6SMSYTWIiVTBY62hKBe4BZIrIJK0k0cT6YK9vkYSzQw3HJ7HrgIHDaGHMM65LXFqfBcU+MBn4FfgC2u7D9V1gLAm0BJju+m4SVHBY6fq6fgX84tp8NPCUiv4lIgzz29yhwr+N7dwGPFeFnUQFKy6or5QYRKQNkOi4FtQM+NMa0LOx7BexvOfCkMcbjKawiUsEYkywikViTCjo4xjtsURw/k/Jtel1SKffUA+aISBCQhnV5rCiOAzNEZJTx/F6OhSJSGQgDXrI5afwEXAqk2xWD8j7tcSillHKLjnEopZRyiyYOpZRSbtHEoZRSyi2aOJRSSrlFE4dSSim3aOJQSinllv8HXeec75xcw0YAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "plt.loglog(sun_obs.wavelength.magnitude, sun_obs.flux.magnitude, label='observed')\n", "plt.loglog(sun_th.wavelength.magnitude, sun_th.flux.magnitude, label='theoretical')\n", "plt.legend();\n", "plt.xlabel('wavelength [{0:s}]'.format(str(sun_obs.wavelength.units)));\n", "plt.ylabel('flux [{0:s}]'.format(str(sun_obs.flux.units)));" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "One can see the differences between the two flavors.\n", "The theoretical spectrum is scaled to match the observed spectrum from 1.5 - 2.5 microns, and then it is used where the observed spectrum ends. The theoretical model of the Sun from Kurucz‘93 atlas using the following parameters when the Sun is at 1 au.\n", "\n", "|log_Z | T_eff | log_g | V$_{Johnson}$ |\n", "|------|-------|-------|---------------|\n", "|+0.0 | 5777 |+4.44 | -26.76 |\n", "\n", "The Sun is also know to have a Johnson V (vega-)magnitude of -26.76 mag.\n", "\n", "Let's verify this." ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "observed 1.84614e+02 -26.7648\n", "theoretical 1.84436e+02 -26.7637\n", "th. 10pc 4.33506e-11 +4.8084\n" ] } ], "source": [ "f = lib['GROUND_JOHNSON_V']\n", "for name, sun in zip(('observed', 'theoretical', 'th. 10pc'), (sun_obs,sun_th, sun_th_10pc)):\n", " flux = f.get_flux(sun.wavelength, sun.flux)\n", " vegamag = f.Vega_zero_mag\n", " print('{0:12s} {1:0.5e} {2:+3.4f}'.format(name, flux.magnitude, -2.5 * np.log10(flux.magnitude) - vegamag))" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": [ "filter_names = ['GROUND_JOHNSON_B', 'GROUND_JOHNSON_V', 'GROUND_BESSELL_J', 'GROUND_BESSELL_K']\n", "filter_names += lib.find('GaiaDR2')\n", "\n", "filters = lib.load_filters(filter_names, lamb=sun_th.wavelength)" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ " GROUND_JOHNSON_B -26.0709 mag\n", " GROUND_JOHNSON_V -26.7637 mag\n", " GROUND_BESSELL_J -27.8970 mag\n", " GROUND_BESSELL_K -28.2706 mag\n", " GaiaDR2_BP -26.5833 mag\n", " GaiaDR2_G -26.9171 mag\n", " GaiaDR2_RP -27.4022 mag\n", " GaiaDR2_weiler_BPbright -26.5835 mag\n", " GaiaDR2_weiler_BPfaint -26.5560 mag\n", " GaiaDR2_weiler_G -26.9115 mag\n", " GaiaDR2_weiler_RP -27.4055 mag\n", " GaiaDR2v2_BP -26.5868 mag\n", " GaiaDR2v2_G -26.9101 mag\n", " GaiaDR2v2_RP -27.3971 mag\n" ] } ], "source": [ "mags = {}\n", "for name, fn in zip(filter_names, filters):\n", " flux = fn.get_flux(sun_th.wavelength, sun_th.flux)\n", " vegamag = fn.Vega_zero_mag\n", " mag = -2.5 * np.log10(flux.magnitude) - vegamag\n", " mags[name] = mag\n", " print('{0:>25s} {1:+3.4f} mag'.format(name, mag))" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ " GROUND_JOHNSON_B - GROUND_JOHNSON_V = 0.6928 mag\n", " GROUND_JOHNSON_V - GROUND_BESSELL_K = 1.5069 mag\n", " GROUND_BESSELL_J - GROUND_BESSELL_K = 0.3735 mag\n", " GaiaDR2_BP - GaiaDR2_RP = 0.8188 mag\n", " GaiaDR2_BP - GaiaDR2_G = 0.3337 mag\n", " GaiaDR2_G - GaiaDR2_RP = 0.4851 mag\n", " GaiaDR2v2_BP - GaiaDR2v2_RP = 0.8103 mag\n", " GaiaDR2v2_BP - GaiaDR2v2_G = 0.3233 mag\n", " GaiaDR2v2_G - GaiaDR2v2_RP = 0.4870 mag\n", " GaiaDR2_weiler_BPbright - GaiaDR2_weiler_RP = 0.8220 mag\n", " GaiaDR2_weiler_BPfaint - GaiaDR2_weiler_RP = 0.8495 mag\n", " GaiaDR2_weiler_BPbright - GaiaDR2_weiler_G = 0.3280 mag\n", " GaiaDR2_weiler_BPfaint - GaiaDR2_weiler_G = 0.3555 mag\n", " GaiaDR2_weiler_G - GaiaDR2_weiler_RP = 0.4940 mag\n" ] } ], "source": [ "colors = (('GROUND_JOHNSON_B', 'GROUND_JOHNSON_V'),\n", " ('GROUND_JOHNSON_V', 'GROUND_BESSELL_K'),\n", " ('GROUND_BESSELL_J', 'GROUND_BESSELL_K'),\n", " ('GaiaDR2_BP', 'GaiaDR2_RP'),\n", " ('GaiaDR2_BP', 'GaiaDR2_G'),\n", " ('GaiaDR2_G', 'GaiaDR2_RP'),\n", " ('GaiaDR2v2_BP', 'GaiaDR2v2_RP'),\n", " ('GaiaDR2v2_BP', 'GaiaDR2v2_G'),\n", " ('GaiaDR2v2_G', 'GaiaDR2v2_RP'),\n", " ('GaiaDR2_weiler_BPbright', 'GaiaDR2_weiler_RP'),\n", " ('GaiaDR2_weiler_BPfaint', 'GaiaDR2_weiler_RP'),\n", " ('GaiaDR2_weiler_BPbright', 'GaiaDR2_weiler_G'),\n", " ('GaiaDR2_weiler_BPfaint', 'GaiaDR2_weiler_G'),\n", " ('GaiaDR2_weiler_G', 'GaiaDR2_weiler_RP'))\n", "\n", "color_values = {}\n", "\n", "for color in colors:\n", " color_values[color] = mags[color[0]] - mags[color[1]]\n", " print('{0:>25s} - {1:<25s} = {2:3.4f} mag'.format(color[0], color[1], mags[color[0]] - mags[color[1]]))\n", " " ] }, { "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.6.8" }, "latex_envs": { "LaTeX_envs_menu_present": true, "autoclose": false, "autocomplete": true, "bibliofile": "biblio.bib", "cite_by": "apalike", "current_citInitial": 1, "eqLabelWithNumbers": true, "eqNumInitial": 1, "hotkeys": { "equation": "ctrl-shift-i" }, "labels_anchors": false, "latex_user_defs": false, "report_style_numbering": false, "user_envs_cfg": false }, "toc": { "base_numbering": 1, "nav_menu": {}, "number_sections": true, "sideBar": true, "skip_h1_title": false, "title_cell": "Table of Contents", "title_sidebar": "Contents", "toc_cell": false, "toc_position": {}, "toc_section_display": true, "toc_window_display": false }, "varInspector": { "cols": { "lenName": 16, "lenType": 16, "lenVar": 40 }, "kernels_config": { "python": { "delete_cmd_postfix": "", "delete_cmd_prefix": "del ", "library": "var_list.py", "varRefreshCmd": "print(var_dic_list())" }, "r": { "delete_cmd_postfix": ") ", "delete_cmd_prefix": "rm(", "library": "var_list.r", "varRefreshCmd": "cat(var_dic_list()) " } }, "types_to_exclude": [ "module", "function", "builtin_function_or_method", "instance", "_Feature" ], "window_display": false } }, "nbformat": 4, "nbformat_minor": 2 }