
# coding: utf-8

# In[2]:

import numpy as np
import netCDF4
import sys
import math 
from pandas import DataFrame 
#sys.path.append("/Users/Beniamino/Documents/OneDrive/MPI-M/iPython_Notebooks/Boreal_Area/") # Edited .zshrc profile with PYTHONPATH 

#from collections import Counter


# In[ ]:




# In[3]:

# execute this only if you want to do the basemap plot

nc005ea = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/\
NEA_50_shift_Bor_MOD44B_V5__VCFFplusVCF__2010__0.05deg__UHAM-ICDC.nc')
ncea = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/\
NEA_50_shift_Bor_MOD44B_V5__VCFFplusVCF__2010__0.5deg__UHAM-ICDC.nc')
nc2ea = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/NEA_50_shift_Bor_MAR_cru_ts3.22.1998.2010.pre.dat.nc')
nc3ea = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/NEA_50_seasonal_soilw.mon.mean.1998-2010.v2.nc')
#nc3ea = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/NEA_shift_min_bor_yearmean_soilw.mon.mean.1998-2010.v2.nc')
nc4ea = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/NEA_50_shift_0_bor_M6H_0.5_tmin.2m.gauss.1998-2010.nc')
nc5ea = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/NEA_50_shift_bor_GLC2000_0.5.nc')
nc6ea = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/NEA_50_shift_0_bor_0.5_PZI.nc')
nc7ea = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/NEA_50_shift_masked_bor_0.5_Year_Fire_Frequency_1996-2012.nc')
# nc8ea = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/NEA_shift_bor_0.5_Burned_1996-2012_SUM.nc')
nc9ea = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/NEA_50_shift_Bor_ReBil_GloElev_05.nc')
nc10ea = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/NEA_50_Boreal_Mean_GDD0_1998-2010_0.5.nc')
nc11ea = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/NEA_50_Bor_Ease25_Northen_Emisphere_Mean_Thawing_Depth_1901-2002.nc')

nc005na = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/\
NA_shift_Bor_MOD44B_V5__VCFFplusVCF__2010__0.05deg__UHAM-ICDC.nc')
ncna = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/\
NA_shift_Bor_MOD44B_V5__VCFFplusVCF__2010__0.5deg__UHAM-ICDC.nc')
nc2na = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/NA_shift_Bor_MAR_cru_ts3.22.1998.2010.pre.dat.nc')
nc3na = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/NA_seasonal_soilw.mon.mean.1998-2010.v2.nc')
#nc3na = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/NA_shift_min_bor_yearmean_soilw.mon.mean.1998-2010.v2.nc')
nc4na = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/NA_shift_0_bor_M6H_0.5_tmin.2m.gauss.1998-2010.nc')
nc5na = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/NA_shift_bor_GLC2000_0.5.nc')
nc6na = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/NA_shift_0_bor_0.5_PZI.nc')
nc7na = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/NA_shift_masked_bor_0.5_Year_Fire_Frequency_1996-2012.nc')
# nc8na = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/NA_shift_bor_0.5_Burned_1996-2012_SUM.nc')
nc9na = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/NA_shift_Bor_ReBil_GloElev_05.nc')
nc10na = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/NA_Boreal_Mean_GDD0_1998-2010_0.5.nc')
nc11na = netCDF4.Dataset('/Users/Beniamino/Documents/Datasets/Boreal/NA_Bor_Ease25_Northen_Emisphere_Mean_Thawing_Depth_1901-2002.nc')

lon_ea = ncea.variables['lon']#.tolist()
lon_na = ncna.variables['lon']

lat_ea = ncea.variables['lat']
lat_na = ncna.variables['lat']

lon005ea = nc005ea.variables['lon']#.tolist()
lon005na = nc005na.variables['lon']

lat005ea = nc005ea.variables['lat']


vcff_ea = ncea.variables['vcff'][0]
vcff005_ea = nc005ea.variables['vcff'][0]

pcp_ea = nc2ea.variables['pre'][0]
sm_ea = nc3ea.variables['soilw']
tm_ea = nc4ea.variables['tmin'][0]-273.15
glc_ea = nc5ea.variables['var'][:]
pzi_ea = nc6ea.variables['pzi'][:]
ff_ea = nc7ea.variables['firefreq'][:]
# ba_ea = nc8ea.variables['burnedarea'][:]
elev_ea = nc9ea.variables['elev'][:]
gdd0_ea = nc10ea.variables['gdd0'][0]
ptd_ea = nc11ea.variables['thw_depth'][0]


vcff_na = ncna.variables['vcff'][0]
vcff005_na = nc005na.variables['vcff'][0]

pcp_na = nc2na.variables['pre'][0]
sm_na = nc3na.variables['soilw']
tm_na = nc4na.variables['tmin'][0]-273.15
glc_na = nc5na.variables['var'][:]
pzi_na = nc6na.variables['pzi'][:]
ff_na = nc7na.variables['firefreq'][:]
# ba_na = nc8na.variables['burnedarea'][:]
elev_na = nc9na.variables['elev'][:]
gdd0_na = nc10na.variables['gdd0'][0]
ptd_na = nc11na.variables['thw_depth'][0]


picctr = 0


# In[4]:

# create new variable of state accordind to the results of the optimisation
vcf_state_na = np.empty(shape=vcff_na.shape)
for i in range(vcff_na.shape[0]):
    for j in range(vcff_na.shape[1]):
        if vcff_na[i][j]<0: 
            vcf_state_na[i][j] = 0  # null
        elif vcff_na[i][j]<20:
            vcf_state_na[i][j] = 1  # bare, the first peak
        #elif vcff_na[i][j]<20:
        #    vcf_state_na[i][j] = 2  # sparse savanna, the second peak
        elif vcff_na[i][j]<45:
            vcf_state_na[i][j] = 3  # savanna, the third peak
        elif vcff_na[i][j]<=100:
            vcf_state_na[i][j] = 4  # forest, the fourth peak
        else:
            vcf_state_na[i][j] = 5  # null

        


# In[5]:

# no need to execute this egain
# create dataset for R with the raveled arrays of data

## IMPORTANT: changes in numpy -> do np.ma.getdata(np.ravel())

lons,lats = np.ma.getdata(np.meshgrid(lon_na,lat_na))
arlons, arlats = np.ma.getdata(np.ravel(lons)), np.ma.getdata(np.ravel(lats))
arsm_na = np.ma.getdata(np.ravel(sm_na[1]))
arvcff_na = np.ma.getdata(np.ravel(vcff_na))
arpcp_na = np.ma.getdata(np.ravel(pcp_na))
artm_na = np.ma.getdata(np.ravel(tm_na))
arglc_na = np.ma.getdata(np.ravel(glc_na))
arpzi_na = np.ma.getdata(np.ravel(pzi_na))
arff_na = np.ma.getdata(np.ravel(ff_na))
#arba_na = np.ravel(ba_na)
arelev_na = np.ma.getdata(np.ravel(elev_na))
argdd0_na = np.ma.getdata(np.ravel(gdd0_na))
arptd_na = np.ma.getdata(np.ravel(ptd_na))
arstate_na = np.ma.getdata(np.ravel(vcf_state_na))


count = 0
index = []
for i in range(len(arsm_na)): 
    if (arsm_na[i]<0 or arpcp_na[i]>4000 or arvcff_na[i]>100 or arglc_na[i]>14 \
        or arpzi_na[i]<0 or arff_na[i]<0 or arelev_na[i]>2000 or arelev_na[i]<0 or argdd0_na[i]<0 or \
       arptd_na[i]<0):
        index.append(i)
        count+=1
data_sm_na = np.delete(arsm_na,index)
data_vcff_na = np.delete(arvcff_na,index)
data_pcp_na = np.delete(arpcp_na,index)
data_tm_na = np.delete(artm_na, index)
data_pzi_na = np.delete(arpzi_na, index)
data_ff_na = np.delete(arff_na, index)
#data_ba_na = np.delete(arba, index)
data_elev_na = np.delete(arelev_na, index)
data_gdd0_na = np.delete(argdd0_na, index)
data_ptd_na = np.delete(arptd_na, index)
data_state_na = np.delete(arstate_na, index)

data_lons, data_lats = np.delete(arlons,index), np.delete(arlats,index)


dataforr_na = np.transpose(np.array([data_vcff_na, data_pcp_na, data_sm_na, data_tm_na, data_pzi_na,\
                                     data_ff_na, data_gdd0_na, data_ptd_na, data_state_na, data_elev_na, data_lons, data_lats])) 
dataframe_na_all = DataFrame(dataforr_na[:,0:9], columns=["VCFF", "MAR", "MASM", "Tmin", "PZI", "FF", "GDD0", "PTD","State"])



# In[9]:

# save the data to a txt that can be reloaded  
#np.savetxt("boreal_na_database.txt", dataforr_na)


# In[4]:

# dataforr_na = np.loadtxt("boreal_na_database.txt")
# dataframe_na_all = DataFrame(dataforr_na[:,0:8], columns=["VCFF", "MAR", "MASM", "Tmin", "PZI", "FF", "GDD0", "State"])


# In[6]:

# no need to execute it again
arvcff005_na = np.ravel(vcff005_na)

count = 0
index = []
for i in range(len(arvcff005_na)): 
    if (arvcff005_na[i]>100 or arvcff005_na[i]<0):
        index.append(i)
        count+=1

data_vcff005_na = np.delete(arvcff005_na,index)
# np.savetxt("boreal_na_vcff_005.txt", data_vcff005_na)


# In[7]:

# data_vcff005_na = np.loadtxt("boreal_na_vcff_005.txt")


# In[8]:

arcsin_vcff005_na = np.arcsin(data_vcff005_na/100.)
arcsin_vcff_na = np.arcsin(data_vcff_na/100.)

#Slice for slicing dataforr skipping State but inlcuding ST (soil texture)
slc = (0,1,2,3,4,5,6,7,8)


# divide datasets according to different state
bare_na = []
sparse_savanna_na = []
savanna_na = []
forest_na = []
null_na = []
for i in range(len(dataforr_na[:,-4])): 
    if (dataforr_na[i,-4] == 0):
        null_na.append(dataforr_na[i,0:8])
    elif (dataforr_na[i,-4] == 1):
        bare_na.append(dataforr_na[i,0:8])
    elif (dataforr_na[i,-4] == 2):
        sparse_savanna_na.append(dataforr_na[i,0:8])
    elif (dataforr_na[i,-4] == 3):
        savanna_na.append(dataforr_na[i,0:8])
    elif (dataforr_na[i,-4] == 4):
        forest_na.append(dataforr_na[i,0:8])
    else:
        null_na.append(dataforr_na[i,0:8])

bare_na = np.array(bare_na)
sparse_savanna_na = np.array(sparse_savanna_na)
savanna_na = np.array(savanna_na)
forest_na = np.array(forest_na)

path_to_databases = "/Users/Beniamino/Documents/OneDrive/MPI-M/iPython_Notebooks/Boreal_Area/"
# In[ ]:
dataforr_na_west = np.loadtxt(path_to_databases+"boreal_na_west_database.txt")
# divide datasets according to different state
bare_na_west = []
sparse_savanna_na_west = []
savanna_na_west = []
forest_na_west = []
null_na_west = []
bare_na_west_coord = []
sava_na_west_coord = []
fore_na_west_coord = []

for i in range(len(dataforr_na_west[:,-4])): 
    if (dataforr_na_west[i,-4] == 0):
        null_na_west.append(dataforr_na_west[i,0:8])
    elif (dataforr_na_west[i,-4] == 1):
        bare_na_west.append(dataforr_na_west[i,0:8])
        bare_na_west_coord.append(dataforr_na_west[i,-2:])
    elif (dataforr_na_west[i,-4] == 2):
        sparse_savanna_na_west.append(dataforr_na_west[i,0:8])
    elif (dataforr_na_west[i,-4] == 3):
        savanna_na_west.append(dataforr_na_west[i,0:8])
        sava_na_west_coord.append(dataforr_na_west[i,-2:])
    elif (dataforr_na_west[i,-4] == 4):
        forest_na_west.append(dataforr_na_west[i,0:8])
        fore_na_west_coord.append(dataforr_na_west[i,-2:])
    else:
        null_na_west.append(dataforr_na_west[i,0:8])


        
bare_na_west = DataFrame(bare_na_west)
bare_na_west[8] = np.array(bare_na_west_coord)[:,-2]
bare_na_west[9] = np.array(bare_na_west_coord)[:,-1]

savanna_na_west = DataFrame(savanna_na_west)
savanna_na_west[8] = np.array(sava_na_west_coord)[:,-2]
savanna_na_west[9] = np.array(sava_na_west_coord)[:,-1]

forest_na_west = DataFrame(forest_na_west)
forest_na_west[8] = np.array(fore_na_west_coord)[:,-2]
forest_na_west[9] = np.array(fore_na_west_coord)[:,-1]
        
bare_na_west = np.array(bare_na_west)
sparse_savanna_na_west = np.array(sparse_savanna_na_west)
savanna_na_west = np.array(savanna_na_west)
forest_na_west = np.array(forest_na_west)



dataforr_na_east = np.loadtxt(path_to_databases+"boreal_na_east_database.txt")#("boreal_na_east_database.txt")
# divide datasets according to different state
bare_na_east = []
sparse_savanna_na_east = []
savanna_na_east = []
forest_na_east = []
null_na_east = []
bare_na_east_coord = []
sava_na_east_coord = []
fore_na_east_coord = []

for i in range(len(dataforr_na_east[:,-4])): 
    if (dataforr_na_east[i,-4] == 0):
        null_na_east.append(dataforr_na_east[i,0:8])
    elif (dataforr_na_east[i,-4] == 1):
        bare_na_east.append(dataforr_na_east[i,0:8])
        bare_na_east_coord.append(dataforr_na_east[i,-2:])
    elif (dataforr_na_east[i,-4] == 2):
        sparse_savanna_na_east.append(dataforr_na_east[i,0:8])
    elif (dataforr_na_east[i,-4] == 3):
        savanna_na_east.append(dataforr_na_east[i,0:8])
        sava_na_east_coord.append(dataforr_na_east[i,-2:])
    elif (dataforr_na_east[i,-4] == 4):
        forest_na_east.append(dataforr_na_east[i,0:8])
        fore_na_east_coord.append(dataforr_na_east[i,-2:])
    else:
        null_na_east.append(dataforr_na_east[i,0:8])
        
bare_na_east = DataFrame(bare_na_east)
bare_na_east[8] = np.array(bare_na_east_coord)[:,-2]
bare_na_east[9] = np.array(bare_na_east_coord)[:,-1]

savanna_na_east = DataFrame(savanna_na_east)
savanna_na_east[8] = np.array(sava_na_east_coord)[:,-2]
savanna_na_east[9] = np.array(sava_na_east_coord)[:,-1]

forest_na_east = DataFrame(forest_na_east)
forest_na_east[8] = np.array(fore_na_east_coord)[:,-2]
forest_na_east[9] = np.array(fore_na_east_coord)[:,-1]


bare_na_east = np.array(bare_na_east)
sparse_savanna_na_east = np.array(sparse_savanna_na_east)
savanna_na_east = np.array(savanna_na_east)
forest_na_east = np.array(forest_na_east)


# In[ ]:
dataforr_na_west_update = np.loadtxt(path_to_databases+"boreal_na_west_database_canadianupdate.txt")
# divide datasets according to different state
bare_na_west_update = []
sparse_savanna_na_west_update = []
savanna_na_west_update = []
forest_na_west_update = []
null_na_west_update = []
bare_na_west_update_coord = []
sava_na_west_update_coord = []
fore_na_west_update_coord = []

for i in range(len(dataforr_na_west_update[:,-4])): 
    if (dataforr_na_west_update[i,-4] == 0):
        null_na_west_update.append(dataforr_na_west_update[i,slc])
    elif (dataforr_na_west_update[i,-4] == 1):
        bare_na_west_update.append(dataforr_na_west_update[i,slc])
        bare_na_west_update_coord.append(dataforr_na_west_update[i,-2:])
    elif (dataforr_na_west_update[i,-4] == 2):
        sparse_savanna_na_west_update.append(dataforr_na_west_update[i,slc])
    elif (dataforr_na_west_update[i,-4] == 3):
        savanna_na_west_update.append(dataforr_na_west_update[i,slc])
        sava_na_west_update_coord.append(dataforr_na_west_update[i,-2:])
    elif (dataforr_na_west_update[i,-4] == 4):
        forest_na_west_update.append(dataforr_na_west_update[i,slc])
        fore_na_west_update_coord.append(dataforr_na_west_update[i,-2:])
    else:
        null_na_west_update.append(dataforr_na_west_update[i,slc])


        
bare_na_west_update = DataFrame(bare_na_west_update)
bare_na_west_update[9] = np.array(bare_na_west_update_coord)[:,-2]
bare_na_west_update[10] = np.array(bare_na_west_update_coord)[:,-1]

savanna_na_west_update = DataFrame(savanna_na_west_update)
savanna_na_west_update[9] = np.array(sava_na_west_update_coord)[:,-2]
savanna_na_west_update[10] = np.array(sava_na_west_update_coord)[:,-1]

forest_na_west_update = DataFrame(forest_na_west_update)
forest_na_west_update[9] = np.array(fore_na_west_update_coord)[:,-2]
forest_na_west_update[10] = np.array(fore_na_west_update_coord)[:,-1]
        
bare_na_west_update = np.array(bare_na_west_update)
sparse_savanna_na_west_update = np.array(sparse_savanna_na_west_update)
savanna_na_west_update = np.array(savanna_na_west_update)
forest_na_west_update = np.array(forest_na_west_update)



dataforr_na_east_update = np.loadtxt(path_to_databases+"boreal_na_east_database_canadianupdate.txt")#("boreal_na_east_database.txt")
# divide datasets according to different state
bare_na_east_update = []
sparse_savanna_na_east_update = []
savanna_na_east_update = []
forest_na_east_update = []
null_na_east_update = []
bare_na_east_update_coord = []
sava_na_east_update_coord = []
fore_na_east_update_coord = []

for i in range(len(dataforr_na_east_update[:,-4])): 
    if (dataforr_na_east_update[i,-4] == 0):
        null_na_east_update.append(dataforr_na_east_update[i,slc])
    elif (dataforr_na_east_update[i,-4] == 1):
        bare_na_east_update.append(dataforr_na_east_update[i,slc])
        bare_na_east_update_coord.append(dataforr_na_east_update[i,-2:])
    elif (dataforr_na_east_update[i,-4] == 2):
        sparse_savanna_na_east_update.append(dataforr_na_east_update[i,slc])
    elif (dataforr_na_east_update[i,-4] == 3):
        savanna_na_east_update.append(dataforr_na_east_update[i,slc])
        sava_na_east_update_coord.append(dataforr_na_east_update[i,-2:])
    elif (dataforr_na_east_update[i,-4] == 4):
        forest_na_east_update.append(dataforr_na_east_update[i,slc])
        fore_na_east_update_coord.append(dataforr_na_east_update[i,-2:])
    else:
        null_na_east_update.append(dataforr_na_east_update[i,slc])
# using a slice like (0,1,2,3,4,5,6,7,9,onetolast,last) the next steps should not be necessary       
bare_na_east_update = DataFrame(bare_na_east_update)
bare_na_east_update[9] = np.array(bare_na_east_update_coord)[:,-2]
bare_na_east_update[10] = np.array(bare_na_east_update_coord)[:,-1]

savanna_na_east_update = DataFrame(savanna_na_east_update)
savanna_na_east_update[9] = np.array(sava_na_east_update_coord)[:,-2]
savanna_na_east_update[10] = np.array(sava_na_east_update_coord)[:,-1]

forest_na_east_update = DataFrame(forest_na_east_update)
forest_na_east_update[9] = np.array(fore_na_east_update_coord)[:,-2]
forest_na_east_update[10] = np.array(fore_na_east_update_coord)[:,-1]


bare_na_east_update = np.array(bare_na_east_update)
sparse_savanna_na_east_update = np.array(sparse_savanna_na_east_update)
savanna_na_east_update = np.array(savanna_na_east_update)
forest_na_east_update = np.array(forest_na_east_update)
        
        
dataforr_na_east_of = np.loadtxt(path_to_databases+"boreal_na_east_onlyfire_database.txt")#("boreal_na_east_database.txt")
# divide datasets according to different state
bare_na_east_of = []
sparse_savanna_na_east_of = []
savanna_na_east_of = []
forest_na_east_of = []
null_na_east_of = []
bare_na_east_of_coord = []
sava_na_east_of_coord = []
fore_na_east_of_coord = []

for i in range(len(dataforr_na_east_of[:,-4])): 
    if (dataforr_na_east_of[i,-4] == 0):
        null_na_east_of.append(dataforr_na_east_of[i,0:8])
    elif (dataforr_na_east_of[i,-4] == 1):
        bare_na_east_of.append(dataforr_na_east_of[i,0:8])
        bare_na_east_of_coord.append(dataforr_na_east_of[i,-2:])
    elif (dataforr_na_east_of[i,-4] == 2):
        sparse_savanna_na_east_of.append(dataforr_na_east_of[i,0:8])
    elif (dataforr_na_east_of[i,-4] == 3):
        savanna_na_east_of.append(dataforr_na_east_of[i,0:8])
        sava_na_east_of_coord.append(dataforr_na_east_of[i,-2:])
    elif (dataforr_na_east_of[i,-4] == 4):
        forest_na_east_of.append(dataforr_na_east_of[i,0:8])
        fore_na_east_of_coord.append(dataforr_na_east_of[i,-2:])
    else:
        null_na_east_of.append(dataforr_na_east_of[i,0:8])


        
bare_na_east_of = DataFrame(bare_na_east_of)
bare_na_east_of[8] = np.array(bare_na_east_of_coord)[:,-2]
bare_na_east_of[9] = np.array(bare_na_east_of_coord)[:,-1]

savanna_na_east_of = DataFrame(savanna_na_east_of)
savanna_na_east_of[8] = np.array(sava_na_east_of_coord)[:,-2]
savanna_na_east_of[9] = np.array(sava_na_east_of_coord)[:,-1]

forest_na_east_of = DataFrame(forest_na_east_of)
forest_na_east_of[8] = np.array(fore_na_east_of_coord)[:,-2]
forest_na_east_of[9] = np.array(fore_na_east_of_coord)[:,-1]


bare_na_east_of = np.array(bare_na_east_of)
sparse_savanna_na_east_of = np.array(sparse_savanna_na_east_of)
savanna_na_east_of = np.array(savanna_na_east_of)
forest_na_east_of = np.array(forest_na_east_of)




dataforr_na_east_nf = np.loadtxt(path_to_databases+"boreal_na_east_nofire_database.txt")#("boreal_na_east_database.txt")
# divide datasets according to different state
bare_na_east_nf = []
sparse_savanna_na_east_nf = []
savanna_na_east_nf = []
forest_na_east_nf = []
null_na_east_nf = []
bare_na_east_nf_coord = []
sava_na_east_nf_coord = []
fore_na_east_nf_coord = []

for i in range(len(dataforr_na_east_nf[:,-4])): 
    if (dataforr_na_east_nf[i,-4] == 0):
        null_na_east_nf.append(dataforr_na_east_nf[i,0:8])
    elif (dataforr_na_east_nf[i,-4] == 1):
        bare_na_east_nf.append(dataforr_na_east_nf[i,0:8])
        bare_na_east_nf_coord.append(dataforr_na_east_nf[i,-2:])
    elif (dataforr_na_east_nf[i,-4] == 2):
        sparse_savanna_na_east_nf.append(dataforr_na_east_nf[i,0:8])
    elif (dataforr_na_east_nf[i,-4] == 3):
        savanna_na_east_nf.append(dataforr_na_east_nf[i,0:8])
        sava_na_east_nf_coord.append(dataforr_na_east_nf[i,-2:])
    elif (dataforr_na_east_nf[i,-4] == 4):
        forest_na_east_nf.append(dataforr_na_east_nf[i,0:8])
        fore_na_east_nf_coord.append(dataforr_na_east_nf[i,-2:])
    else:
        null_na_east_nf.append(dataforr_na_east_nf[i,0:8])


        
bare_na_east_nf = DataFrame(bare_na_east_nf)
bare_na_east_nf[8] = np.array(bare_na_east_nf_coord)[:,-2]
bare_na_east_nf[9] = np.array(bare_na_east_nf_coord)[:,-1]

savanna_na_east_nf = DataFrame(savanna_na_east_nf)
savanna_na_east_nf[8] = np.array(sava_na_east_nf_coord)[:,-2]
savanna_na_east_nf[9] = np.array(sava_na_east_nf_coord)[:,-1]

forest_na_east_nf = DataFrame(forest_na_east_nf)
forest_na_east_nf[8] = np.array(fore_na_east_nf_coord)[:,-2]
forest_na_east_nf[9] = np.array(fore_na_east_nf_coord)[:,-1]


bare_na_east_nf = np.array(bare_na_east_nf)
sparse_savanna_na_east_nf = np.array(sparse_savanna_na_east_nf)
savanna_na_east_nf = np.array(savanna_na_east_nf)
forest_na_east_nf = np.array(forest_na_east_nf)



dataforr_ea = np.loadtxt(path_to_databases+"boreal_ea_database.txt")#("boreal_na_east_database.txt")
# divide datasets according to different state
bare_ea = []
sparse_savanna_ea = []
savanna_ea = []
forest_ea = []
null_ea = []
bare_ea_coord = []
sava_ea_coord = []
fore_ea_coord = []

for i in range(len(dataforr_ea[:,-4])): 
    if (dataforr_ea[i,-4] == 0):
        null_ea.append(dataforr_ea[i,0:8])
    elif (dataforr_ea[i,-4] == 1):
        bare_ea.append(dataforr_ea[i,0:8])
        bare_ea_coord.append(dataforr_ea[i,-2:])
    elif (dataforr_ea[i,-4] == 2):
        sparse_savanna_ea.append(dataforr_ea[i,0:8])
    elif (dataforr_ea[i,-4] == 3):
        savanna_ea.append(dataforr_ea[i,0:8])
        sava_ea_coord.append(dataforr_ea[i,-2:])
    elif (dataforr_ea[i,-4] == 4):
        forest_ea.append(dataforr_ea[i,0:8])
        fore_ea_coord.append(dataforr_ea[i,-2:])
    else:
        null_ea.append(dataforr_ea[i,0:8])


        
bare_ea = DataFrame(bare_ea)
bare_ea[8] = np.array(bare_ea_coord)[:,-2]
bare_ea[9] = np.array(bare_ea_coord)[:,-1]

savanna_ea = DataFrame(savanna_ea)
savanna_ea[8] = np.array(sava_ea_coord)[:,-2]
savanna_ea[9] = np.array(sava_ea_coord)[:,-1]

forest_ea = DataFrame(forest_ea)
forest_ea[8] = np.array(fore_ea_coord)[:,-2]
forest_ea[9] = np.array(fore_ea_coord)[:,-1]


bare_ea = np.array(bare_ea)
sparse_savanna_ea = np.array(sparse_savanna_ea)
savanna_ea = np.array(savanna_ea)
forest_ea = np.array(forest_ea)


# In[ ]:
dataforr_ea_west = np.loadtxt(path_to_databases+"boreal_ea_west_database.txt")
# divide datasets according to different state
bare_ea_west = []
sparse_savanna_ea_west = []
savanna_ea_west = []
forest_ea_west = []
null_ea_west = []
bare_ea_west_coord = []
sava_ea_west_coord = []
fore_ea_west_coord = []

for i in range(len(dataforr_ea_west[:,-4])): 
    if (dataforr_ea_west[i,-4] == 0):
        null_ea_west.append(dataforr_ea_west[i,0:8])
    elif (dataforr_ea_west[i,-4] == 1):
        bare_ea_west.append(dataforr_ea_west[i,0:8])
        bare_ea_west_coord.append(dataforr_ea_west[i,-2:])
    elif (dataforr_ea_west[i,-4] == 2):
        sparse_savanna_ea_west.append(dataforr_ea_west[i,0:8])
    elif (dataforr_ea_west[i,-4] == 3):
        savanna_ea_west.append(dataforr_ea_west[i,0:8])
        sava_ea_west_coord.append(dataforr_ea_west[i,-2:])
    elif (dataforr_ea_west[i,-4] == 4):
        forest_ea_west.append(dataforr_ea_west[i,0:8])
        fore_ea_west_coord.append(dataforr_ea_west[i,-2:])
    else:
        null_ea_west.append(dataforr_ea_west[i,0:8])


        
bare_ea_west = DataFrame(bare_ea_west)
bare_ea_west[8] = np.array(bare_ea_west_coord)[:,-2]
bare_ea_west[9] = np.array(bare_ea_west_coord)[:,-1]

savanna_ea_west = DataFrame(savanna_ea_west)
savanna_ea_west[8] = np.array(sava_ea_west_coord)[:,-2]
savanna_ea_west[9] = np.array(sava_ea_west_coord)[:,-1]

forest_ea_west = DataFrame(forest_ea_west)
forest_ea_west[8] = np.array(fore_ea_west_coord)[:,-2]
forest_ea_west[9] = np.array(fore_ea_west_coord)[:,-1]
        
bare_ea_west = np.array(bare_ea_west)
sparse_savanna_ea_west = np.array(sparse_savanna_ea_west)
savanna_ea_west = np.array(savanna_ea_west)
forest_ea_west = np.array(forest_ea_west)



dataforr_ea_east = np.loadtxt(path_to_databases+"boreal_ea_east_database.txt")#("boreal_ea_east_database.txt")
# divide datasets according to different state
bare_ea_east = []
sparse_savanna_ea_east = []
savanna_ea_east = []
forest_ea_east = []
null_ea_east = []
bare_ea_east_coord = []
sava_ea_east_coord = []
fore_ea_east_coord = []

for i in range(len(dataforr_ea_east[:,-4])): 
    if (dataforr_ea_east[i,-4] == 0):
        null_ea_east.append(dataforr_ea_east[i,slc])
    elif (dataforr_ea_east[i,-4] == 1):
        bare_ea_east.append(dataforr_ea_east[i,slc])
        bare_ea_east_coord.append(dataforr_ea_east[i,-2:])
    elif (dataforr_ea_east[i,-4] == 2):
        sparse_savanna_ea_east.append(dataforr_ea_east[i,slc])
    elif (dataforr_ea_east[i,-4] == 3):
        savanna_ea_east.append(dataforr_ea_east[i,slc])
        sava_ea_east_coord.append(dataforr_ea_east[i,-2:])
    elif (dataforr_ea_east[i,-4] == 4):
        forest_ea_east.append(dataforr_ea_east[i,slc])
        fore_ea_east_coord.append(dataforr_ea_east[i,-2:])
    else:
        null_ea_east.append(dataforr_ea_east[i,slc])
        
bare_ea_east = DataFrame(bare_ea_east)
bare_ea_east[9] = np.array(bare_ea_east_coord)[:,-2]
bare_ea_east[10] = np.array(bare_ea_east_coord)[:,-1]

savanna_ea_east = DataFrame(savanna_ea_east)
savanna_ea_east[9] = np.array(sava_ea_east_coord)[:,-2]
savanna_ea_east[10] = np.array(sava_ea_east_coord)[:,-1]

forest_ea_east = DataFrame(forest_ea_east)
forest_ea_east[9] = np.array(fore_ea_east_coord)[:,-2]
forest_ea_east[10] = np.array(fore_ea_east_coord)[:,-1]


bare_ea_east = np.array(bare_ea_east)
sparse_savanna_ea_east = np.array(sparse_savanna_ea_east)
savanna_ea_east = np.array(savanna_ea_east)
forest_ea_east = np.array(forest_ea_east)


# In[ ]:
dataforr_ea_west_russia = np.loadtxt(path_to_databases+"boreal_ea_west_russia_database.txt")
# divide datasets according to different state
bare_ea_west_russia = []
sparse_savanna_ea_west_russia = []
savanna_ea_west_russia = []
forest_ea_west_russia = []
null_ea_west_russia = []
bare_ea_west_russia_coord = []
sava_ea_west_russia_coord = []
fore_ea_west_russia_coord = []

for i in range(len(dataforr_ea_west_russia[:,-4])): 
    if (dataforr_ea_west_russia[i,-4] == 0):
        null_ea_west_russia.append(dataforr_ea_west_russia[i,slc])
    elif (dataforr_ea_west_russia[i,-4] == 1):
        bare_ea_west_russia.append(dataforr_ea_west_russia[i,slc])
        bare_ea_west_russia_coord.append(dataforr_ea_west_russia[i,-2:])
    elif (dataforr_ea_west_russia[i,-4] == 2):
        sparse_savanna_ea_west_russia.append(dataforr_ea_west_russia[i,slc])
    elif (dataforr_ea_west_russia[i,-4] == 3):
        savanna_ea_west_russia.append(dataforr_ea_west_russia[i,slc])
        sava_ea_west_russia_coord.append(dataforr_ea_west_russia[i,-2:])
    elif (dataforr_ea_west_russia[i,-4] == 4):
        forest_ea_west_russia.append(dataforr_ea_west_russia[i,slc])
        fore_ea_west_russia_coord.append(dataforr_ea_west_russia[i,-2:])
    else:
        null_ea_west_russia.append(dataforr_ea_west_russia[i,slc])


        
bare_ea_west_russia = DataFrame(bare_ea_west_russia)
bare_ea_west_russia[9] = np.array(bare_ea_west_russia_coord)[:,-2]
bare_ea_west_russia[10] = np.array(bare_ea_west_russia_coord)[:,-1]

savanna_ea_west_russia = DataFrame(savanna_ea_west_russia)
savanna_ea_west_russia[9] = np.array(sava_ea_west_russia_coord)[:,-2]
savanna_ea_west_russia[10] = np.array(sava_ea_west_russia_coord)[:,-1]

forest_ea_west_russia = DataFrame(forest_ea_west_russia)
forest_ea_west_russia[9] = np.array(fore_ea_west_russia_coord)[:,-2]
forest_ea_west_russia[10] = np.array(fore_ea_west_russia_coord)[:,-1]
        
bare_ea_west_russia = np.array(bare_ea_west_russia)
sparse_savanna_ea_west_russia = np.array(sparse_savanna_ea_west_russia)
savanna_ea_west_russia = np.array(savanna_ea_west_russia)
forest_ea_west_russia = np.array(forest_ea_west_russia)



dataforr_ea_europe = np.loadtxt(path_to_databases+"boreal_ea_europe_database.txt")#("boreal_ea_europe_database.txt")
# divide datasets according to different state
bare_ea_europe = []
sparse_savanna_ea_europe = []
savanna_ea_europe = []
forest_ea_europe = []
null_ea_europe = []
bare_ea_europe_coord = []
sava_ea_europe_coord = []
fore_ea_europe_coord = []

for i in range(len(dataforr_ea_europe[:,-4])): 
    if (dataforr_ea_europe[i,-4] == 0):
        null_ea_europe.append(dataforr_ea_europe[i,slc])
    elif (dataforr_ea_europe[i,-4] == 1):
        bare_ea_europe.append(dataforr_ea_europe[i,slc])
        bare_ea_europe_coord.append(dataforr_ea_europe[i,-2:])
    elif (dataforr_ea_europe[i,-4] == 2):
        sparse_savanna_ea_europe.append(dataforr_ea_europe[i,slc])
    elif (dataforr_ea_europe[i,-4] == 3):
        savanna_ea_europe.append(dataforr_ea_europe[i,slc])
        sava_ea_europe_coord.append(dataforr_ea_europe[i,-2:])
    elif (dataforr_ea_europe[i,-4] == 4):
        forest_ea_europe.append(dataforr_ea_europe[i,slc])
        fore_ea_europe_coord.append(dataforr_ea_europe[i,-2:])
    else:
        null_ea_europe.append(dataforr_ea_europe[i,slc])
        
bare_ea_europe = DataFrame(bare_ea_europe)
bare_ea_europe[9] = np.array(bare_ea_europe_coord)[:,-2]
bare_ea_europe[10] = np.array(bare_ea_europe_coord)[:,-1]

savanna_ea_europe = DataFrame(savanna_ea_europe)
savanna_ea_europe[9] = np.array(sava_ea_europe_coord)[:,-2]
savanna_ea_europe[10] = np.array(sava_ea_europe_coord)[:,-1]

forest_ea_europe = DataFrame(forest_ea_europe)
forest_ea_europe[9] = np.array(fore_ea_europe_coord)[:,-2]
forest_ea_europe[10] = np.array(fore_ea_europe_coord)[:,-1]


bare_ea_europe = np.array(bare_ea_europe)
sparse_savanna_ea_europe = np.array(sparse_savanna_ea_europe)
savanna_ea_europe = np.array(savanna_ea_europe)
forest_ea_europe = np.array(forest_ea_europe)