{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# TESS 2015 - Introduction to Solar Data Analysis in Python" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "# SunPy!\n", "## Author: Steven Christe\n", "Email: steven.christe@nasa.gov" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Through tutorials and presentations we will demonstrate how free and open-source Python and SunPy library can be used to analyze solar data. Depending on interest, dinner may be ordered after the main presentation (roughly an hour) and a hands-on help session will take place for the remainder of the evening. Installing the following software is recommended but not required: Anaconda Python and SunPy." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Schedule: 19:00 to 22:00\n", "+ 18:00-18:15 pm: Organization (chatroom, organize dinner)?\n", "+ 18:15-19:00 pm: Intro to SunPy (presentation)\n", "+ 19:00-19:15 pm: Break\n", "+ 19:15-20:15: Dinner & Installation Workshop\n", "+ 20:15-22:00: SunPy Workshop" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## What is SunPy?\n", "### A community-developed, free and open-source solar data analysis environment for Python." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "+ website: [http://www.sunpy.org](http://www.sunpy.org)
\n", "+ documentation: [http://docs.sunpy.org](http://docs.sunpy.org)
\n", "+ code (version control!): [https://github.com/sunpy/sunpy](https://github.com/sunpy/sunpy)\n", "+ Mailing list: https://groups.google.com/forum/#!forum/sunpy**\n", "+ IRC: #sunpy on freenode.net [web client](https://kiwiirc.com/client/irc.freenode.net/#SunPy)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "SunPy is built upon foundational libraries which enable scientific computing in Python which includes\n", "\n", "+ [NumPy](http://numpy.org)\n", "+ [SciPy](http://scipy.org)\n", "+ [matplotlib](http://matplotlib.org)\n", "+ [AstroPy](http://astropy.org)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "

Supported observations\n", "

\n", " \n", " \n", " \n", "

Supported data retrieval services \n", "

\n", "

\n", " \n", "

Supported file formats include\n", "

\n", "

" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# What is this webby notebooky magic?!\n", "## ipython notebook (now known as jupyter)\n", "similar to matlab, mathematica, maple" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "Some setup..." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "import numpy as np\n", "import matplotlib.pyplot as plt\n", "import matplotlib\n", "%matplotlib inline\n", "matplotlib.rc('savefig', dpi=120)\n", "import warnings\n", "warnings.simplefilter(\"ignore\", Warning)\n", "from matplotlib import dates" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "import sunpy" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "sunpy.system_info()" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "SunPy version (stable) 0.5" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Let's study a flare!" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Searching for events in the HEK" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "from sunpy.net import hek\n", "client = hek.HEKClient()\n", "tstart, tend = '2014/01/01 00:00:00', '2014/01/02 00:00:00'\n", "result = client.query(hek.attrs.Time(tstart, tend), \n", " hek.attrs.EventType('FL'), \n", " hek.attrs.FRM.Name=='SSW Latest Events')\n", "len(result)" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "result[0]" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "for res in result:\n", " print(res.get('fl_goescls'))" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "result = client.query(hek.attrs.Time(tstart, tend), \n", " hek.attrs.EventType('FL'), \n", " hek.attrs.FRM.Name=='SSW Latest Events', \n", " hek.attrs.FL.GOESCls>'M')" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "len(result)" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "result" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can find out when this event occured" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "result[0].get('event_peaktime')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "and where it occurred" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "result[0].get('hpc_coord')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Lightcurves!\n", "## Let's look at the GOES curve for this event.\n", "First some time manipulation!" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": true }, "outputs": [], "source": [ "from sunpy.time import TimeRange, parse_time\n", "from datetime import timedelta" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": true }, "outputs": [], "source": [ "tmax = parse_time(result[0].get('event_peaktime'))" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "tmax" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": true }, "outputs": [], "source": [ "tr = TimeRange(tmax - timedelta(minutes=30), tmax + timedelta(minutes=30))" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "tr" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "from sunpy.lightcurve import GOESLightCurve\n", "goes = GOESLightCurve.create(tr)" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "goes.peek()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The data is stored in a standard place" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "goes.data" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This is a pandas dataframe! Provides lots of additional functionality. For example" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "print('The max flux is {flux:2.5f} at {time}'.format(flux=goes.data['xrsb'].max(), time=goes.data['xrsb'].idxmax()))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Compares well to the official max from the HEK" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "str(tmax)" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "goes.peek()\n", "plt.axhline(goes.data['xrsb'].max())\n", "plt.axvline(goes.data['xrsb'].idxmax())" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Meta data is also stored in a standard place" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "goes.meta" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This is a dictionary like the hek results so..." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "goes.meta.get('COMMENT')" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "goes.data.resample('10s', how='mean')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Solar Images in SunPy" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "SunPy has a `Map` type that supports 2D images, it makes it simple to read data in from any filetype supported in `sunpy.io` which is currently FITS, JPEG2000 and ANA files. You can also create maps from any `(data, metadata)` pair." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let's download an AIA image of this flare from the vso" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": true }, "outputs": [], "source": [ "from sunpy.net import vso\n", "client=vso.VSOClient()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Then do a search." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "recs = client.query(vso.attrs.Time(tr), vso.attrs.Instrument('AIA'))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "It works, now lets see how many results we have!" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "recs.num_records()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "That's way too many!" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "So that is every image that SDO/AIA had on that day. Let us reduce that amount.\n", "\n", "To do this, we will limit the time search for a specify a wavelength." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": true }, "outputs": [], "source": [ "recs = client.query(vso.attrs.Time('2014/01/01 18:52:08', '2014/01/01 18:52:15'), \n", " vso.attrs.Instrument('AIA'),\n", " vso.attrs.Wave(171,171))" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "recs.num_records()" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "recs.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let's also grab another wavelength for later." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": true }, "outputs": [], "source": [ "recs = client.query(vso.attrs.Time('2014/01/01 18:52:08', '2014/01/01 18:52:15'), \n", " vso.attrs.Instrument('AIA'),\n", " vso.attrs.Wave(94,171))" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "recs.num_records()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let's download this data!" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "f = client.get(recs, methods = ('URL-FILE_Rice')).wait()" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "f" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "For SunPy the top level name-space is kept clean. Importing SunPy does not give you access to much. You need to import specific names. SciPy is the same. \n", "\n", "So, the place to start here is with the core SunPy data object. It is called Map. " ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": true }, "outputs": [], "source": [ "from sunpy.map import Map" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "aia = Map(f[1])" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "aia" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Maps contain both the image data and the metadata associated with the image, this metadata currently does not deviate much from the standard FITS WCS keywords, but presented in a instrument-independent manner." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "aia.peek()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# SunPy Maps!\n", "Maps are the same for each file, **regardless of source**. It does not matter if the source is SDO or SOHO, for example.\n", "\n", "The most used attributes are as follows (some of them will look similar to NumPy's Array):" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "aia.data" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The data (stored in a numpy array)" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "type(aia.data)" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "aia.mean(),aia.max(),aia.min()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Because it is just a numpy array you have access to all of those function" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The standard deviation" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "aia.data.std()" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "aia.data.shape" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The original metadata (stored in a dictionary)" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "aia.meta" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "aia.meta.keys()" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "aia.meta.get('rsun_obs')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We also provide quick access to some key metadata values as object variables (these are shortcuts)" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "print(aia.date, aia.coordinate_system, aia.detector, aia.dsun)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Maps also provide some nice map specific functions such as submaps. Let's zoom in on the flare location which was given to us by the HEK." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "result[0].get('hpc_coord')" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "point = [665.04, -233.4096]\n", "dx = 50\n", "dy = 50\n", "xrange = [point[0] - dx, point[0] + dx]\n", "yrange = [point[1] - dy, point[1] + dy]\n", "aia.submap(xrange,yrange).peek()\n", "plt.plot(point[0], point[1], '+')\n", "plt.xlim(xrange)\n", "plt.ylim(yrange)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The default image scale is definitely not right. Let's fix that so we can see the flare region better." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "smap = aia.submap(xrange,yrange)\n", "\n", "import matplotlib.colors as colors\n", "norm = colors.Normalize(0, 3000)\n", "smap.plot(norm=norm)\n", "plt.plot(point[0], point[1], '+')\n", "smap.draw_grid(grid_spacing=1)\n", "plt.colorbar()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Composite Maps\n", "Let's plot the two channels we downloaded from AIA together! Composite map is the way to do this. Let's check out the other map we got." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "aia131 = Map(f[0])\n", "aia131" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "smap131 = aia131.submap(xrange, yrange)\n", "smap131.peek()" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "norm = colors.Normalize(0, 4000)\n", "smap131.plot(norm=norm)\n", "plt.colorbar()" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": true }, "outputs": [], "source": [ "smap171 = smap" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "compmap = Map(smap171, smap131, composite=True)" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "levels = np.arange(0,100,5)\n", "print(levels)" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "compmap.set_levels(1, levels, percent=True)\n", "compmap.set_mpl_color_normalizer(0, norm)\n", "compmap.set_colors(1, plt.cm.Reds)" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "compmap.plot(norm=norm)\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Some other topics..." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Solar Constants" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "from sunpy.sun import constants as solar_constants" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "solar_constants.mass" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "print(solar_constants.mass)" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "(solar_constants.mass/solar_constants.volume).cgs" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "solar_constants.volume + solar_constants.density" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Not all constants have a short-cut assigned to them (as above). The rest of the constants are stored in a dictionary. The following code grabs the dictionary and gets all of the keys.:" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "solar_constants.physical_constants.keys()" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "type(solar_constants.mass)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "These are astropy constants which are a subclass of Quantities (numbers with units) which are a great idea." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": true }, "outputs": [], "source": [ "from astropy import units as u" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "u.keV" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "u.keV.decompose()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This has been a simple overview of AstroPy units. IF you want to read more, see http://astropy.readthedocs.org/en/latest/units/." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## More References\n", "\n", "+ [Python for Data Analysis](http://shop.oreilly.com/product/0636920023784.do) by Wes McKinney is great\n", "+ [Lectures on scientific computing with Python](http://jrjohansson.github.io) Great lectures on scientific computing in Python" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Consider [contributing](http://sunpy.org/contribute/) to SunPy!\n", "\n", "+ Provide feedback on the [mailing list](https://groups.google.com/forum/#!forum/sunpy)\n", "+ Report [bugs](https://github.com/sunpy/sunpy/issues)\n", "+ Provide Code (see our [developer guide](http://sunpy.readthedocs.org/en/stable/dev.html))" ] } ], "metadata": { "kernelspec": { "display_name": "Python 2", "language": "python", "name": "python2" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 2 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython2", "version": "2.7.9" } }, "nbformat": 4, "nbformat_minor": 0 }