# ---
# jupyter:
#   jupytext:
#     text_representation:
#       extension: .py
#       format_name: light
#       format_version: '1.5'
#       jupytext_version: 1.14.4
#   kernelspec:
#     display_name: Python 3
#     language: python
#     name: python3
# ---

# ## Compile anomalies and NPC position into one file for analysis
# #### created by sol, nov 26 019
# ### New data from chelle - nov 23, 2019 - v6

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import xarray as xr

# +
# npc position
ds = xr.open_dataset('../../biforcation_figures_chelle/biforcation_data165_monthly_aviso2.nc')
ds2=ds.lat.sel(time=slice('1993-01-01','2018-12-31'))

# anomalies
clim=ds2.groupby('time.month').mean('time')
ds3=ds2.groupby('time.month')-clim
npc165 = ds3.values

ds = xr.open_dataset('../../biforcation_figures_chelle/biforcation_data135_monthly_aviso2.nc')
ds2=ds.lat.sel(time=slice('1993-01-01','2018-12-31'))
fecha=ds2.time.values
# anomalies
clim=ds2.groupby('time.month').mean('time')
ds3=ds2.groupby('time.month')-clim
npc135 = ds3.values
# -

plt.plot(fecha,npc165,label='165')
plt.plot(fecha,npc135,label='135')
plt.legend()
plt.title('NPC Position')
plt.show()

# +
boxl=['NPC','Bif','Cal2','Ak1','Ak2','NPCv','Bifv', 'Cal1','Cal3']
# currents
ds = xr.open_dataset('../../timeseries_data/oscar1993-01-01data_minus_clim_v6.nc')
# select V Cal1, V Ak1 for this analysis
ds2=ds.sel(time=slice('1993-01-01','2018-12-31'))
vcak = ds2.v[3,:].values 
vcac = ds2.v[7,:].values
plt.plot(fecha,vcak,label='CAK')
plt.plot(fecha,vcac,label='CAC')
plt.legend()
plt.title('Oscar V currents')
plt.show()
# SST
ds = xr.open_dataset('../../timeseries_data/sst1993-01-01data_minus_clim_v6.nc')
# select V Cal1, V Ak1 for this analysis
ds2=ds.sel(time=slice('1993-01-01','2018-12-31'))
tcak = ds2.analysed_sst[3,:].values
tcac = ds2.analysed_sst[7,:].values
plt.plot(fecha,tcak, label='CAK')
plt.plot(fecha,tcac, label='CAC')


t165 = ds2.analysed_sst[5,:].values
t135 = ds2.analysed_sst[6,:].values
plt.plot(fecha,t165, label='165W')
plt.plot(fecha,t135, label='135W')
plt.title('SST')

plt.legend()
plt.show()
# -

ds

var = ['NPC165', 'NPC135','V_CAC','VCAK']
adt1 = xr.DataArray(npc165, dims=['time'], coords=[fecha])
adt2 = xr.DataArray(npc135, dims=['time'], coords=[fecha])
adt3 = xr.DataArray(vcac, dims=['time'], coords=[fecha])
adt4 = xr.DataArray(vcak, dims=['time'], coords=[fecha])
adt5 = xr.DataArray(tcac, dims=['time'], coords=[fecha])
adt6 = xr.DataArray(tcak, dims=['time'], coords=[fecha])
adt7 = xr.DataArray(t165, dims=['time'], coords=[fecha])
adt8 = xr.DataArray(t135, dims=['time'], coords=[fecha])
ads = xr.Dataset({'NPC165': adt1, 'NPC135':adt2, 'V_CAC': adt3, 'V_CAK':adt4, 'SST_CAC': adt5, 'SST_CAK':adt6, 'SST_NPC': adt7, 'SST_Bif':adt8})
ads.to_netcdf('../data/bifdatams_long_26nov2019.nc',mode='w')

ads

# seasonal anomalies
adseas = ads.resample(time='3M',label='left',closed='left',loffset='1M').mean()
adseas = adseas.sel(time=slice('1993-01-01','2018-12-31'))
adseas.to_netcdf('../data/bifdatams_seas_long_26nov2019.nc',mode='w')

dfseas = adseas.to_dataframe()
dfseas = dfseas.rename(columns={'NPC165':'NPC@165W','NPC135':'NPC@135W','V_CAC':'CalCv','V_CAK':'AkCv','SST_NPC':'SST@165W','SST_Bif':'SST@135W'})
dfseas = dfseas.drop(['SST_CAC','SST_CAK'],axis=1)
dfseas.to_csv('../data/MS_data_feb2020.csv')

plt.figure(figsize=(16,4))
adseas.SST_NPC.plot()
ads.SST_NPC.plot()
plt.grid()

plt.figure(figsize=(16,4))
adseas.SST_Bif.plot()
ads.SST_Bif.plot()
plt.grid()

ds = xr.open_dataset('../../timeseries_data/sst1993-01-01data_minus_clim_v6.nc')
ds

312/3



StarcoderdataJupyterScriptsDedupFiltered