#! /usr/bin/env python
# -*- coding:utf8 -*-
import matplotlib.pyplot as plt
import numpy as np
from scipy import stats
from math import *
plt.style.use('classic')
# download data from ftp://ftp.marine.csiro.au/pub/legresy/gmsl_files/CSIRO_Alt.csv
data = []
with open('CSIRO_Alt.csv', 'r') as csv:
for txt in csv.readlines():
try:
line = [i.strip() for i in txt.split(',')]
if len(line) == 3:
data.append(line)
except Exception:
pass
years0, h0 = [], []
years1, h1 = [], []
for d in data:
try:
y, h = float(d[0]), float(d[1])
years0.append(y); h0.append(h)
except Exception:
pass
try:
y, h = float(d[0]), float(d[2])
years1.append(y); h1.append(h)
except Exception:
pass
years0 = np.array(years0)
slope, intercept, r_value, p_value, std_err = stats.linregress(years0, h0)
plt.rc('font', size=15)
fig = plt.figure(figsize=(800 / 90.0, 540 / 90.0), dpi=72)
plt.plot(years0, h0, '-', color='#00eeee', linewidth=2, label='Monthly')
plt.plot(years1, h1, '-', color='#0000dd', linewidth=3, label='3-month running mean')
plt.plot(years0[[0,-1]], intercept+slope*np.array(years0[[0,-1]]), '-',
color='red', linewidth=1.5, label='Trend = {:.1f} mm/year'.format(slope))
plt.grid(True)
plt.xlim(1992, ceil(max(years0)/1.)*1.)
plt.gca().ticklabel_format(useOffset=False)
plt.ylim(-50, 59)
plt.gca().ticklabel_format(useOffset=False)
plt.xlabel('Year')
plt.ylabel('Global Mean Sea Level (mm)')
plt.title('GMSL from TOPEX/Poseidon, Jason-1, Jason-2 and Jason-3 satellite altimeter data', size=14)
plt.legend(loc='upper left', prop={'size':14}, borderaxespad=1.2)
plt.text(0.47, 0.058, '''Seasonal signal removed
Inverse barometer correction applied
GIA correction applied''', size=14, transform=plt.gca().transAxes)
plt.tight_layout()
plt.savefig('Alt_gmsl_seas_rem.svg')