"hint: checkout the [matplotlib examples](http://matplotlib.org/examples/index.html)"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"# %load files/sol_sinc.py"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The SciPy framework builds on top NumPy and provides a large number of numerical algorithms for working with data. Some of the topics that SciPy covers are:\n",
"\n",
"* Special functions ([scipy.special](http://docs.scipy.org/doc/scipy/reference/special.html))\n",
"In many cases, you will want to import specific names from `scipy` subpackages. However, as a start, it is helpful to do the following import:"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"import scipy as sp"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Approach"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"One of the most important skills in data science is to be able to find Python functions and classes in a module and learn how to use them yourself. Here are some recommended steps on how to go about this:\n",
"\n",
"* Find the online documentation for the package you are using.\n",
"* Try to find the subpackage or even the function that looks like will do the job.\n",
"* Import the module, function or class and use tab completion and `?` to explore it.\n",
"* Try using the function or class for an extremely simple case where you know the answer.\n",
"* Then try using for your real problem."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Data Fitting\n",
"\n",
"Let's randomize a exponential decay function and try to do a linear fit.\n",
"Adapted from [Python for scientific use, Part II: Data analysis](http://linuxgazette.net/115/andreasen.html) \n",
"by Anders Andreasen. See also a more modern approach to the subject using [Pandas](http://pandas.pydata.org/), the Python Data Analysis Library, on [this Josh Hemann notebook](http://nbviewer.ipython.org/gist/jhemann/4569783)\n",
"\n",
"In the next example we will use the Fast Fourier Transform (FFT) in order to transform time-dependent data into the frequency domain. By doing so, it is possible to analyse if any predominant frequencies exists - i.e. if there is any periodicity in the data. \n",
"\n",
"The example data we will use is the sunspot activity measurements from year 1700 to year 2014 provided by the Sunspot index and Long-term Solar Observation. The data set is the yearly mean total sunspot observations available via HTTP.\n",
"\n",
"The datafile is also included in `files/SN_y_tot_V2.0.txt`. The sunspot data have been used for illustrating the power of FFT with respect to finding a periodicity in sunspot activity in various computer languages, e.g. Matlab and BASIC. "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Retrieve data\n",
"\n",
"Solar spot data are available online from the [datafiles section](http://www.sidc.be/silso/datafiles) of the Sunspot index and Long-term Solar Observation webpage. We will use the yearly mean total sunspot number [plain text version.]( versionhttp://www.sidc.be/silso/DATA/SN_y_tot_V2.0.txt)\n",
"\n",
"A quick web search point us to the `urllib` python module and to a very simple retrival procedure."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### FFT\n",
"\n",
"Let's take a look at the [SciPy Fast Fourier Transform module manual](http://docs.scipy.org/doc/scipy/reference/tutorial/fftpack.html), and straight from there...\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"In order to construct a periodogram, i.e. a graph of power vs. frequency, we first compute the power of the FFT signal which is simply the **FFT signal squared**. We only need the part of the signal ranging from zero to a frequency equal to the Nyquist frequency, which is equal to half the maximum frequency, since frequencies above the Nyquist frequency correspond to negative frequencies. The frequency range is calculated from 0-N/2 as N/(2T) where N is the number of samples and T is the sampling time."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Which in years becomes"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Solar spots period"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"###Scipy direct method"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"After doing the exercise, while looking for information, I found that scipy allows us to directly compute the periodgram with [`scipy.signal.periodgram`](http://scipy.github.io/devdocs/generated/scipy.signal.periodogram.html)"
