i1 mag trends+spline

# -*- coding: UTF-8 -*-

#Modifed 2012/11/16

import numpy as np
import matplotlib as mpl
import matplotlib.pyplot as plt
from numpy import *
import sys
from scipy.interpolate import spline

plt.clf()

ffile = sys.argv[1]

data = genfromtxt(ffile, unpack=True)

N   = data[0][:]
mod = data[1][:]
fi  = data[2][:]
i1  = data[3][:]

X = N

Y_01 = mod
Y_02 = i1

#Magnitude linear trend
coefficients_mod = polyfit(X, Y_01, 1)
polynomial_mod   = poly1d(coefficients_mod)
xs_mod           = arange(0, max(N), 100)
ys_mod           = polynomial_mod(xs_mod)

#Magnitude spline 4
spline_coefficients_mod = polyfit(X, Y_01, 4)
spline_polynomial_mod   = poly1d(spline_coefficients_mod)
spline_xs_mod           = arange(0, max(N), 100)
spline_ys_mod           = spline_polynomial_mod(spline_xs_mod)

#i1 linear trend
i1_lt_coefficients_i1  = polyfit(X, Y_02, 1)
i1_lt_polynomial_i1    = poly1d(i1_lt_coefficients_i1)
i1_lt_xs_i1            = arange(0, max(N), 100)
i1_lt_ys_i1            = i1_lt_polynomial_i1(i1_lt_xs_i1)

#i1 spline 4
i1_sp_coefficients_i1  = polyfit(X, Y_02, 4)
i1_sp_polynomial_i1    = poly1d(i1_sp_coefficients_i1)
i1_sp_xs_i1            = arange(0, max(N), 100)
i1_sp_ys_i1            = i1_sp_polynomial_i1(i1_sp_xs_i1)


mpl.rcParams['figure.figsize'] = (8.0, 6.0)

#Magnitude
line_mod,         = plt.plot(X, Y_01, 'ro-', markevery = (0, 5), linewidth=0.3, markersize = 0, label ='Magnitude')
#Linear trend of magnitude
line_spline_mod   = plt.plot(xs_mod, ys_mod,  'r--',  linewidth=1.0, label ='Linear trend of magnitude')
#Spline of magnitude
line_spline_mod   = plt.plot(spline_xs_mod, spline_ys_mod,  'r-',  linewidth=3.0, label ='4th spline of magnitude')

#I1
line_i1           = plt.plot(X, Y_02, 'bo-', markevery = (0, 5), linewidth=0.3, markersize = 0, label = 'Total variability(I1)')
#Linear trend I1
line_lt_i1        = plt.plot(i1_lt_xs_i1, i1_lt_ys_i1,  'b--',  linewidth=1.0, label ='Linear trend of I1')
#Spline I1
line_lt_i1        = plt.plot(i1_sp_xs_i1, i1_sp_ys_i1,  'b-',  linewidth=3.0, label ='4th spline of I1')

plt.title(u'Magnitude and variability (1979-2006)')

plt.legend(loc = 'best')

ax_01 = plt.axes()
ax_01.grid(color = 'black')
ax_01.set_xlabel(u'[years]')
ax_01.set_ylabel(u'[cm/s]')

plt.axis((min(N)-300, max(N)+400, min(mod), max(i1)))

labels = [item.get_text() for item in ax_01.get_xticklabels()]

labels[1] = '1979'
labels[2] = '1984'
labels[3] = '1989'
labels[4] = '1994'
labels[5] = '1999'
labels[6] = '2007'

ax_01.set_xticklabels(labels)

plt.savefig(ffile[0]+'_I1_mod.png', dpi=300, format = 'png')

print 'Correlation (magnitude,i1): %2.1f' % (corrcoef(mod, i1)[0,1])