#!/usr/bin/env python
# coding: utf-8

# In[1]:


# Play with Beaching #


# In[164]:


import datetime as dt
import matplotlib.pyplot as plt
import numpy as np
from pathlib import Path
import xarray as xr

from salishsea_tools import viz_tools


# In[44]:


mesh = xr.open_dataset('/home/sallen/MEOPAR/grid/mesh_mask201702.nc')
land = mesh.tmask[0, 0]
mesh


# In[206]:


# Kits
imin, imax = 448, 460
jmin, jmax = 321, 345
fig, ax = plt.subplots(1, 1)
land[imin:imax, jmin:jmax].plot(ax=ax)
viz_tools.set_aspect(ax);
mymask = np.zeros_like(land)
mymask[imin:imax, jmin:jmax] = 1


# In[210]:


# Orcas
imin, imax = 320, 327
jmin, jmax = 284, 308
fig, ax = plt.subplots(1, 1)
land[imin:imax, jmin:jmax].plot(ax=ax)
viz_tools.set_aspect(ax);
mymask = np.zeros_like(land)
mymask[imin:imax, jmin:jmax] = 1


# In[208]:


# Lummi
imin, imax = 305, 320
jmin, jmax = 337, 355
fig, ax = plt.subplots(1, 1)
land[imin:imax, jmin:jmax].plot(ax=ax)
viz_tools.set_aspect(ax);
mymask = np.zeros_like(land)
mymask[imin:imax, jmin:jmax] = 1


# values = np.ma.masked_array(np.ones_like(land[imin:imax, jmin:jmax]), mask=1-land[imin:imax, jmin:jmax])

# np.ma.notmasked_edges(values)

# In[193]:


minoil = 5


# In[212]:


colours = {'bunker': 'blue',
           'other': 'blue',
           'akns': 'navy',
           'diesel': 'skyblue',
           'gas': 'skyblue',
           'jet': 'skyblue',
           'dilbit': 'navy',
           }
colours['bunker']
data.OilType.item()


# In[176]:


month_colours = ['slateblue', 'slateblue', 'lawngreen', 'lawngreen', 'indianred', 'indianred',
                 'goldenrod', 'goldenrod', 'goldenrod', 'darkorchid', 'darkorchid', 'slateblue']


# In[214]:


scale = 5
scount = 0
ncount = 0
fig, axs = plt.subplots(1, 2, figsize=(15, 7))
for ax in axs:
    ax.pcolormesh(mesh.nav_lon, mesh.nav_lat, land[:-1,:-1], cmap='BuGn_r')
direct = Path('/data/sallen/results/MIDOSS/beaching_files/')
for file in direct.glob('*'):
    data = xr.open_dataset(file)
#    month = data.Spilldatetime.values.astype('datetime64[M]').astype(int) % 12
    if (np.ma.masked_array(data.Beaching_Volume, mask=data.Beaching_Volume < minoil/1000.) * 
        mymask[1:-1, 1:-1]).sum(axis=0).sum(axis=0) > 0 :
        
        axs[1].scatter(data.SpillLon, data.SpillLat, marker='o', c=colours[data.OilType.item()], 
                   s=np.log(1000*data.SpillVolume.item())*scale,
               edgecolors='k', zorder=3)
        scount = scount + 1
    else:
        axs[0].scatter(data.SpillLon, data.SpillLat, marker='x', c=colours[data.OilType.item()], 
                   s=np.log(1000*data.SpillVolume.item())*scale)
        ncount = ncount + 1
        
    data.close()
for ax in axs:
    ax.set_xlim(-123.6, -122.35)
    ax.set_ylim(48.35, 49.4)
    ax.plot([mesh.nav_lon[imin, jmin], mesh.nav_lon[imin, jmax], mesh.nav_lon[imax, jmax], 
             mesh.nav_lon[imax, jmin], mesh.nav_lon[imin, jmin]],
            [mesh.nav_lat[imin, jmin], mesh.nav_lat[imin, jmax], mesh.nav_lat[imax, jmax], 
             mesh.nav_lat[imax, jmin], mesh.nav_lat[imin, jmin]],
            color='orange')
    viz_tools.set_aspect(ax, coords='map');
axs[0].set_title(f'{ncount} Spills that do *not* Result in Oil on Coastline in Orange Box')
axs[1].set_title(f'{scount} Spills that do Result in Oil on Coastline in Orange Box');


# In[37]:


data


# In[171]:


mydate = data.Spilldatetime.values
mydate


# In[168]:


dt.datetime.fromtimestamp(mydate)


# In[175]:


data.Spilldatetime.values.astype('datetime64[M]').astype(int) % 12 + 1


# In[173]:


months


# In[ ]: