"""
Plots Arctic daily sea ice extent from June 2002-present using JAXA metadata

Website   : https://ads.nipr.ac.jp/vishop/vishop-extent.html
Author    : Zachary M. Labe
Date      : 15 May 2016
"""

### Import modules
import numpy as np
import matplotlib.pyplot as plt
import matplotlib
import datetime
import urllib.request
import urllib as UL

### Directory and time
directoryfigure = './Figures/'
now = datetime.datetime.now()
currentmn = str(now.month)
currentdy = str(now.day-1)
currentyr = str(now.year)
currenttime = currentmn + '_' + currentdy + '_' + currentyr

### Load url
url = 'https://ads.nipr.ac.jp/vishop.ver1/data/graph/plot_extent_n_v2.csv'

### Read file
raw_data = urllib.request.urlopen(url)
dataset = np.genfromtxt(raw_data, skip_header=0,delimiter=",",)

### Set missing data to nan
dataset[np.where(dataset==-9999)] = np.nan

### Variables
month     = dataset[1:,0]        # 1-12, nan as month[0]
day       = dataset[1:,1]        # 1-31, nan as day[0]
mean1980  = dataset[1:,2]        # km^2, nan as mean1980[0]
mean1990  = dataset[1:,3]        # km^2, nan as mean1990[0]
mean2000  = dataset[1:,4]        # km^2, nan as mean2000[0]
years     = dataset[1:,5:]

doy       = np.arange(0,len(day),1)

### Change units to million km^2
years = years/1e6

### Recent day of current year
currentyear = years[:,-1] 
lastday = now.timetuple().tm_yday -1
currentice = currentyear[lastday]
currentanom = currentice - (mean1980[lastday]/1e6)

# Leap year
currentyear[59] = currentyear[58]

### Changes in the last day and week
weekchange = currentice - currentyear[lastday-7]
daychange = currentice - currentyear[lastday-1]

###############################################################################
###############################################################################
###############################################################################
### Plot figure
matplotlib.rc('savefig', facecolor='black')
matplotlib.rc('axes', edgecolor='white')
matplotlib.rc('xtick', color='white')
matplotlib.rc('ytick', color='white')
matplotlib.rc('axes', labelcolor='white')
matplotlib.rc('axes', facecolor='black')
plt.rc('text',usetex=True)
plt.rc('font',**{'family':'sans-serif','sans-serif':['Avant Garde']})

fig = plt.figure()
ax = plt.subplot(111) 

### Adjust axes in time series plots 
def adjust_spines(ax, spines):
    for loc, spine in ax.spines.items():
        if loc in spines:
            spine.set_position(('outward', 0))
        else:
            spine.set_color('none')  
    if 'left' in spines:
        ax.yaxis.set_ticks_position('left')
    else:
        ax.yaxis.set_ticks([])

    if 'bottom' in spines:
        ax.xaxis.set_ticks_position('bottom')
    else:
        ax.xaxis.set_ticks([])

oldaverage = currentyear.copy()
oldaverage[lastday:] = currentyear[lastday]

### 2000s average
average2000s = mean2000.copy()
average2000s[lastday:] = mean2000[lastday]
average2000s = average2000s/1e6
oldmin = np.where(mean2000 == np.min(mean2000))[0]

### 1990s average
average1990s = mean1990.copy()
average1990s[lastday:] = mean1990[lastday]
average1990s = average1990s/1e6

### 1980s average
average1980s = mean1980.copy()
average1980s[lastday:] = mean1980[lastday]
average1980s = average1980s/1e6

difference = (oldmin - lastday)[0]

### Are we below decadal climatological min?
if (currentyear[lastday]*1e6) < np.nanmin(mean1980):
    print( True, '1980')
if (currentyear[lastday]*1e6) < np.nanmin(mean1990):
    print(True, '1990')
if (currentyear[lastday]*1e6) < np.nanmin(mean2000):
    print(True, '2000')
    
### Calculate record low SIE
recordlow = np.empty((years.shape[0]))
for i in range(years.shape[0]):
    if years[i,-1] == np.nanmin(years[i,:]):
        recordlow[i] = 1.
    else: 
        recordlow[i] = 0.

### Begin plot
plt.plot(doy,years[:,:],color='w',linewidth=0.15,
         linestyle='-',alpha=0.7)

bar = ax.plot(doy,currentyear,linewidth=2.9,zorder=3,
              color='darkorange',)
                  
plt.scatter(doy[lastday],currentyear[lastday],
            s=20,color='darkorange',zorder=4)
plt.scatter(doy[lastday],mean2000[lastday]/1e6,
            s=20,color='dodgerblue',zorder=11)   
plt.scatter(doy[lastday],mean1990[lastday]/1e6,
            s=20,color='c',zorder=11)  
plt.scatter(doy[lastday],mean1980[lastday]/1e6,
            s=20,color='darkmagenta',zorder=11)             