分组再平均

分组再平均

# -*- coding: utf-8 -*-

"""
Created on  Dec 17 15:30:00 2017

@author: grug in xi'an
"""
import os
import sys
import matplotlib.pyplot as plt
import pandas  as pd
import numpy  as np
import matplotlib.dates as dates
import prettyplotlib as ppl
from prettyplotlib import mpl
from prettyplotlib import brewer2mpl
def search(s, path=os.path.abspath('.')):  #os.path.abspath(path)：绝对路径

    for z in os.listdir(path):
        if os.path.isdir(path + os.path.sep + z):  # os.path.sep:路径分隔符 linux下就用这个了’/’
            #print('Currnet:', path)
            path2 = os.path.join(path, z) #；os.path.join(): 常用来链接路径
            #print('future:', path2)
            search(s, path2)
        elif os.path.isfile(path + os.path.sep + z): #检验给出的路径是否是一个文件：os.path.isfile()来自 <http://blog.csdn.net/devil_2009/article/details/7941241>
            if s in z:
                #print(os.path.join(path, z))
                filenames.append(os.path.join(path, z))
                #with open(path + os.path.sep + z, 'r') as fr:
                #    with open('save.txt', 'a') as fw:
                #        fw.write(path + '\t' + fr.read())
def plotData(ax,x,y,ptype, pcolor,plabel):
    #for chinese!
    #plt.rcParams['font.sans-serif']=['msyh']
    #plt.rcParams['axes.unicode_minus'] = False
    #配置一下坐标刻度等
    #xlabel=np.linspace(0,1,25)
    # You can specify a rotation for the tick labels in degrees or with keywords.
    times=pd.to_datetime(x,format='%Y-%m-%d_%H:%M:%S')
    #ax.plot_date(times.to_pydatetime(), y, ptype,label=plabel,alpha=0.5)
    ax.xaxis.set_minor_locator(dates.DayLocator())
    ax.xaxis.set_minor_formatter(dates.DateFormatter('%d'))
    ax.xaxis.set_major_locator(dates.MonthLocator())
    ax.xaxis.set_major_formatter(dates.DateFormatter('%b'))
    #ax.yaxis.set_ticks_position('left')
    #ax.xaxis.set_ticks_position('bottom')

    ax.xaxis.grid(True, which="minor",color='lightgray', linestyle='-.', linewidth=0.2)
    ax.yaxis.grid(True, which="major",color='lightgray', linestyle='-.', linewidth=0.2)


    #plt.xticks(rotation='30')
    #plt.xticks(fontsize = 5)
    #plt.yticks(fontsize = 5)
    #We change the fontsize of minor ticks label
    plt.tick_params(axis='both', which='major', labelsize=10)
    plt.tick_params(axis='both', which='minor', labelsize=8)#just for the last tick
    # Pad margins so that markers don't get clipped by the axes
    #plt.margins(0.2)
    #plt.xlim(X.min() * 1.1, X.max() * 1.1)
    #plt.ylim(C.min() * 1.1, C.max() * 1.1)
    #plt.axes(fontsize=5)
    #plt.legend()  # 让图例生效
    ppl.plot(ax,times.to_pydatetime(), y,ptype,color=pcolor,label=plabel,alpha=0.9,markersize=1.5,linewidth=1.2) # 加上label参数添加图例
    #plt.tight_layout()
    #pl.show()#让绘制的图像在屏幕上显示出来；注释后才能叠加到一起
    #配置一下坐标刻度等
def getmax(row):
    vmax=mean(row[''],)
    return vmax

value = sys.argv[1]
chem= '../chem/nanjing.csv'
chem100= '../100_chem/nanjing.csv'
chem10= '../10_chem/nanjing.csv'
chem02= '../02_chem/nanjing.csv'
obsdata= './Pre_nanjing_use.csv'
#obsdata= 'observation/pm_yrd_use.csv'

data=pd.read_csv(chem)#数量不定的空格符可以用正则表达式\s+表示
data100=pd.read_csv(chem100)#数量不定的空格符可以用正则表达式\s+表示
data10=pd.read_csv(chem10)#数量不定的空格符可以用正则表达式\s+表示
data02=pd.read_csv(chem02)#数量不定的空格符可以用正则表达式\s+表示
obs=pd.read_csv(obsdata)#数量不定的空格符可以用正则表达式\s+表示
data.set_index(['BJC','STA'], inplace=True)
print(data)
data100.set_index(['BJC','STA'], inplace=True)
print(data100)
data02.set_index(['BJC','STA'], inplace=True)
print(data02)
data10.set_index(['BJC','STA'], inplace=True)
print(data10)
obs.set_index(['BJC','Stid'], inplace=True)
data=data.iloc[data.index.isin(['JN', 'PK', 'LH','LS','GC'], level=1)]
data02=data02.iloc[data02.index.isin(['JN', 'PK', 'LH','LS','GC'], level=1)]
data10=data10.iloc[data10.index.isin(['JN', 'PK', 'LH','LS','GC'], level=1)]
data100=data100.iloc[data100.index.isin(['JN', 'PK', 'LH','LS','GC'], level=1)]
obs=obs.iloc[obs.index.isin(['JN', 'PK', 'LH','LS','GC'], level=1)]
#station=data.index.levels[1]
station= [ u'JN',u'PK',u'LH',
           u'LS',u'GC'
           ]
print(station)
data100=data100['PRE'].groupby(level=[0]).apply(lambda x: x.mean())
print(data100.index)
print(data100)
data10=data10['PRE'].groupby(level=[0]).apply(lambda x: x.mean())
data02=data02['PRE'].groupby(level=[0]).apply(lambda x: x.mean())
data=data['PRE'].groupby(level=[0]).apply(lambda x: x.mean())
obs=obs['PRE'].groupby(level=[0]).apply(lambda x: x.mean())
fig, ax= plt.subplots(1, sharex=True,figsize=(8,6))#sharey=True
#fig.subplots_adjust(left=0.1, bottom=0.08, right=0.9, top=0.96,hspace=0.2, wspace=0.3)
#axes=(ax1, ax2,ax3,ax4,ax5)
plotData(ax,data.index,data,'*-','#e78ac3',u'100%')#注意下站点对不对应
plotData(ax,data02.index,data02,'*-','#8da0cb',u'50%')#注意下站点对不对应
plotData(ax,data10.index,data10,'*-',"#fc8d62",u'10%')#注意下站点对不对应
plotData(ax,data100.index,data100,'*-',"#66c2a5",u'1%')#注意下站点对不对应
plotData(ax,obs.index,obs,'*-','k',u'obs')#注意下站点对不对应
ax.set_title("NanJing",fontsize=8,loc="left")
plt.legend(bbox_to_anchor=(0.8,0.8), loc="lower left",
          borderaxespad=0, ncol=1)
if value[0:2]=='CC':
    fig.text(0.02, 0.5, r'CCN($cm^{-3}$)', va='center', rotation='vertical')
elif value[0:2]=='pm':
    fig.text(0.02, 0.5, r'$PM_{2.5}(\mu g/m^3 $)', va='center', rotation='vertical')
elif value[0:2]=='PR':
    fig.text(0.02, 0.5, r'Precipitation($mm/h$)', va='center', rotation='vertical',fontsize=10)
else:
    print(value)
fig.text(0.5, 0.02, 'Time(BJT)', ha='center',fontsize=10)
#ax.xaxis.set_ticks_position('bottom')
plt.savefig(value+'_5bplot_d04.png',dpi=200,bbox_inches='tight')