# -*- coding: utf-8 -*-
"""
Created on Wed Oct 18 16:35:40 2017
读取白洁师姐土壤水实测数据
@author: zhang
"""
import pandas as pd
import matplotlib.pyplot as plt
import scipy.stats as stats
import time
datafilepath = 'G:/XJ_drought/measured_data/BJ_201405-201409.xlsx'
#rdmifilepath = 'G:/XJ_drought/Try_1st/Verification_BJ_point/RDMI_BJ_2014_1st.xlsx'
rdmifilepath = 'G:/XJ_drought/Try_2nd/RDMI_BJ_2014_2nd.xlsx'
measure_data = pd.read_excel(datafilepath,na_values=['0','99.9'])
rdmi_data = pd.read_excel(rdmifilepath)
rdmi_data['Date'] = rdmi_data['Date'].dt.date
rdmi_data.rename(columns={'Date':'time'},inplace=True)
rdmi_data['RDMIpercent'] = rdmi_data['value']*100
select_time = datetime.time(10,0)
data_10 = measure_data[map(lambda x:x.time()==select_time,measure_data['time'])]
data_10_20cm = data_10['20cm'].fillna(data_10['20cm'].mean())
data_10_noNan = data_10.dropna(axis=0)
#data_20, = plt.plot(data_10_noNan['time'],data_10_noNan['20cm'])
#data_40, = plt.plot(data_10_noNan['time'],data_10_noNan['40cm'],linestyle='--',marker='^')
#data_60, = plt.plot(data_10_noNan['time'],data_10_noNan['60cm'],linestyle='-.',marker='s')
#data_80, = plt.plot(data_10_noNan['time'],data_10_noNan['80cm'],linestyle=':',marker='+')
#plt.legend([data_20,data_40,data_60,data_80],[u'20cm',u'40cm',u'60cm',u'80cm'],fontsize=12,loc='upper left')
#plt.xlabel('Time',fontsize=16)
#plt.ylabel('Water content(%)',fontsize=16)
#plt.xticks(fontsize = 12,rotation = 17)
#plt.yticks(fontsize = 12)
#data_10.to_excel('G:/XJ_drought/measured_data/data_BJ_10.xlsx')
data_10_20cm_noNan = data_10[['time','20cm']].dropna(axis=0)
data_10_20cm_noNan['time'] = data_10_20cm_noNan['time'].dt.date
#plt.plot(data_10_20cm_noNan['time'],data_10_20cm_noNan['20cm'])
data_10_20cm_noNan['severe'] = ((data_10_20cm_noNan['20cm'].max()-data_10_20cm_noNan['20cm'])/
(data_10_20cm_noNan['20cm'].max()-data_10_20cm_noNan['20cm'].min()))*100
contract = pd.merge(data_10_20cm_noNan,rdmi_data,how='inner')
#contract.to_excel('G:/XJ_drought/Try_2nd/Contract_BJ_point_2014.xlsx')
contract['RDMIpercent'].corr(contract['severe'])
contract['RDMIpercent'].cov(contract['20cm'])
r, p=stats.pearsonr(contract['RDMIpercent'],contract['20cm'])
print(r,p)
plt.scatter(contract['RDMIpercent'],contract['severe'])
plt.xlim((0,100))
plt.ylim((0,100))
