import numpy as np
import netCDF4
import datetime as dt
import pandas as pd
import pratic as p


# period data
year=2012
yvec=str(year)


# zones data
zones=['Northern Europe (Dfc)','Central Europe (Cfb)','Eastern Europe (Dfb)','Mediterranean Contour (Csa)', 'Arabia (BWh)', 'South China (Cwa)', 'South Brasil (Cfa)', 'Central America (Aw)','Boreal North America (Dfc)','USA (Bsk,Cfa)','Boreal Asian steppe (Dfb,BSk)','Tropical Asia (Af)','Australia (Bwh,Cfa)','North Africa (BWh)','Equatorial Africa (Aw,Am)','South Africa (Bsh,Bsk)','South America (Cfc,Csa,BSk)','West Asia (Cwb,Cfa,Dwa,Bwk)','Turkey (Csa,Csb)','New Zealand (Cfb)','Mexico (Bsh,Aw)','Ukraine (Dfb)']
ind_zones=[31,32,33,34,37,42,51,50,0,30,40,43,44,46,47,48,52,65,63,58,64]
nb_zones=len(ind_zones) #21


# Input files
obs_d = "/home/satellites7/bdelorme/ESACCI-SM/daily/res05/ESACCI-SOILMOISTURE-L3S-SSMV-COMBINED-%(year)s.nc"
obs_m = "/home/satellites7/bdelorme/ESACCI-SM/monthly/res05/ESACCI-SOILMOISTURE-L3S-SSMV-COMBINED-%(year)s.nc"
mod_d = "/home/satellites7/bdelorme/ORCHIDEE/daily/res05/4c/e2o_cnrs_wrr1_glob30_day_SoilMoist_%(year)s.nc"
mod_m = "/home/satellites7/bdelorme/ORCHIDEE/monthly/res05/4c/e2o_cnrs_wrr1_glob30_day_SoilMoist_%(year)s.nc"
mask = "/home/satellites7/bdelorme/masks/philippe/regions_mask_extend_land_mask05.nc"


# graph data
dates = []
ndays = [31, 29 if year%4==0 else 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31]
for month in range(12):
    for day in range(ndays[month]):
        dates.append(dt.datetime(year=year, month=month+1, day=day+1))


# get data
print "obs daily"
nc = netCDF4.Dataset(obs_d % dict(year=year),"r")
map_obs_d=nc.variables["sm"][:]
nc.close()
print "mod daily"
nc = netCDF4.Dataset(mod_d % dict(year=year),"r")
map_mod_d=nc.variables["sm"][:]
nc.close()
print "obs monthly"
nc = netCDF4.Dataset(obs_m % dict(year=year),"r")
map_obs_m=nc.variables["sm"][:]
nc.close()
print "mod monthly"
nc = netCDF4.Dataset(mod_m % dict(year=year),"r")
map_mod_m=nc.variables["sm"][:]
nc.close()
map_mod_d[map_mod_d>1]=np.nan
map_obs_d[map_obs_d>1]=np.nan
map_mod_m[map_mod_m>1]=np.nan
map_obs_m[map_obs_m>1]=np.nan
map_mod_d=map_mod_d+map_obs_d-map_obs_d
map_obs_d=map_obs_d+map_mod_d-map_mod_d
map_mod_m=map_mod_m+map_obs_m-map_obs_m
map_obs_m=map_obs_m+map_mod_m-map_mod_m


# Spatial autocorrelation
days=[15,45,75,105,135,165,195,225,255,285,315,345]
nlag=30
lmax=3000
print "semivario"
#p.semivariogram(MAT,nlag,lmax)
p.semivario(map_obs_d,map_mod_d,map_obs_m,map_mod_m,mask,nb_zones,zones,yvec,'FD',ind_zones,days,nlag,lmax)
