image.png

image.png

image.png

交通流實(shí)驗(yàn)如上

車輛到達(dá)分布

1027 17 2.33207715978 符合負(fù)二項(xiàng)分布 .png

散點(diǎn)圖

1031 10 11[6, 42, 1297.8333333333333, 367.0, 0.88408265213442316].png

解析出的數(shù)據(jù)

image.png

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
sns.set_style('darkgrid')
#計(jì)算高峰小時(shí)
global end,start,vds
A=[]
B=[]
C=[]
D=[]
M=[]
S=[]
def init(start,end,vds):#初始化 數(shù)據(jù)錄入部分
    df = {}
    didi = []
    for i in range(0, end - start):
        df[i] = pd.read_excel(
            'C:\\Users\Administrator\Desktop\Whitemud Drive地感線圈數(shù)據(jù)（2015年8月）\\' + str(vds) + '-08-' + str(
                start + i) + '.xls', sheetname='Sheet3')
        didi.append(len(df[i]))
    length = int(min(didi) / 100) * 100
    for i in range(0, end - start):
        df[i] = df[i][0:length]
        df[i] = df[i].fillna(0)
    global K,Q,V
    K = np.zeros((length, end - start))
    Q = np.zeros((length, end - start))
    V = np.zeros((length, end - start))
    for i in range(0, end - start):
        for j in range(0, length):
            if df[i]['K 2'][j] != 0 and df[i]['K 3'][j] != 0 and df[i]['K 1'][j] != 0:
                K[j, i] = (df[i]['K 1'][j] + df[i]['K 2'][j] + df[i]['K 3'][j]) / 3
                Q[j, i] = (df[i]['V 1'][j] + df[i]['V 2'][j] + df[i]['V 3'][j]) / 3
                V[j, i] = (df[i]['S 1'][j] + df[i]['S 2'][j] + df[i]['S 3'][j]) / 3
            if (df[i]['K 2'][j] != 0 and df[i]['K 3'][j] == 0 and df[i]['K 1'][j] == 0) or (
                        df[i]['K 2'][j] == 0 and df[i]['K 3'][j] != 0 and df[i]['K 1'][j] == 0) or (
                        df[i]['K 2'][j] == 0 and df[i]['K 3'][j] == 0 and df[i]['K 1'][j] != 0):
                K[j, i] = ((df[i]['K 1'][j] + df[i]['K 2'][j] + df[i]['K 3'][j]))
                Q[j, i] = ((df[i]['V 1'][j] + df[i]['V 2'][j] + df[i]['V 3'][j]))
                V[j, i] = ((df[i]['S 1'][j] + df[i]['S 2'][j] + df[i]['S 3'][j]))
            if (df[i]['K 2'][j] != 0 and df[i]['K 3'][j] != 0 and df[i]['K 1'][j] == 0) or (
                        df[i]['K 2'][j] == 0 and df[i]['K 3'][j] != 0 and df[i]['K 1'][j] != 0) or (
                        df[i]['K 2'][j] != 0 and df[i]['K 3'][j] == 0 and df[i]['K 1'][j] != 0):
                K[j, i] = (df[i]['K 1'][j] + df[i]['K 2'][j] + df[i]['K 3'][j]) / 2
                Q[j, i] = (df[i]['V 1'][j] + df[i]['V 2'][j] + df[i]['V 3'][j]) / 2
                V[j, i] = (df[i]['S 1'][j] + df[i]['S 2'][j] + df[i]['S 3'][j]) / 2
    K[K > 80] = 80



def PHT():#計(jì)算高峰小時(shí)
    Q1=Q[:,0]/180
    for i in range(len(Q1),0,-1):
        if i%3==0 and i%6==0 and i%15==0:
            break
    Q1=Q1[0:i]
    Z=np.array(Q1).reshape((int(i/3),3))
    sum=Z.sum(axis=1)

    if sum.sum()<=10:
        return 0
    maxhour=0
    max15=0
    for i in range(0,len(Z)-60):
        emm=Z[i:i + 60].sum()
        if emm>maxhour:
            maxhour=emm
            starttime=i
    for i in range(0,len(Z)-15):
        emm=Z[i:i + 15].sum()
        if emm>max15:
            max15=emm

    # print(starttime,maxhour,max15)
    # print(maxhour/(4*max15))
    # print(int(starttime/60),(starttime%60))

    # temp=[int(starttime/60),starttime%60,maxhour,max15,maxhour/(max15*4)]
    # time=str(int(starttime/60))+':'+str(starttime%60)
    # C.append(time)
    # A.append(vds)
    # B.append(start)
    # return temp

    return starttime




def modal(min):#檢驗(yàn)分布函數(shù)
    various=min.var()/min.mean()
    if various>1.2:
        print('符合負(fù)二項(xiàng)分布')
        return 1,various
    if various<1.2 and various>1 :
        print('符合泊松分布')
        return 2,various
    if  various < 1:
        print('二項(xiàng)分布')
        return 3,various




def arrival():#車輛到達(dá)分布
    #計(jì)算最大公約數(shù)

    Q1 = Q[:, 0] / 180
    for i in range(len(Q1),0,-1):
        if i%3==0 and i%6==0 and i%15==0:
            break
    Q1=Q1[0:i]
    #計(jì)算頻次
    min1=np.array(Q1[:,0]).reshape((int(i/3),3))
    min2 = np.array(Q1[:,0]).reshape((int(i/6), 6))
    min5 = np.array(Q1[:,0]).reshape((int(i/15) , 15))
    min1=min1.sum(axis=1)
    min2 = min2.sum(axis=1)
    min5 = min5.sum(axis=1)



    if min5.sum()<10:
        return 0

    M.append(min1.mean())
    S.append(min1.var())
    min1=min1[min1>5]
    min2 = min2[min2 >5]
    min5 = min5[min5 >5]
    #驗(yàn)證分布
    temp,various=modal(min1)
    modal(min2)
    modal(min5)
    #畫圖并保存
    plt.subplot(311)
    sns.distplot(min1, color='#ff8000')
    plt.subplot(312)
    sns.distplot(min2, color='#ff8000')
    plt.subplot(313)
    sns.distplot(min5, color='#ff8000')

    #
    # if temp==1:
    #     plt.savefig('E:\\arrival\負(fù)二項(xiàng)分布\\'+str(vds)+' '+str(start)+'  '+str(various)+' 符合負(fù)二項(xiàng)分布'+'  '+'.png')
    # if temp==2:
    #     plt.savefig('E:\\arrival\泊松分布\\'+str(vds)+' '+str(start)+'  '+str(various)+' 符合泊松分布'+'  '+'.png')
    # if temp==3:
    #     plt.savefig('E:\\arrival\二項(xiàng)分布\\'+str(vds)+' '+str(start)+'  '+str(various)+' 二項(xiàng)分布'+'  '+'.png')
    #plt.close()

    D.append(temp)




def hourarrival():#高峰小時(shí)車輛到達(dá)分布
    starttime=PHT()
    #print(starttime*3)
    endtime=starttime*3+180
    Q1=Q[starttime*3:endtime]/180
    #print(Q1[:,0])
    min1 = np.array(Q1[:, 0]).reshape((int(len(Q1)/ 3), 3))
    min2 = np.array(Q1[:, 0]).reshape((int(len(Q1)/ 6), 6))
    min5 = np.array(Q1[:, 0]).reshape((int(len(Q1)/ 15), 15))
    min1 = min1.sum(axis=1)
    min2 = min2.sum(axis=1)
    min5 = min5.sum(axis=1)

    if min5.sum() < 10:
        return 0

    M.append(min1.mean())
    S.append(min1.var())
    min1 = min1[min1 > 1]
    min2 = min2[min2 > 1]
    min5 = min5[min5 > 1]
    # 驗(yàn)證分布
    #temp, various = modal(min1)
    #modal(min2)
    temp,various=modal(min5)
    # 畫圖并保存
    # plt.subplot(311)
    # sns.distplot(min1, color='#ff8000')
    # plt.subplot(312)
    # sns.distplot(min2, color='#ff8000')
    # plt.subplot(313)
    # sns.distplot(min5, color='#ff8000')
    # #plt.show()
    # if temp==1:
    #     plt.savefig('E:\\arrival\hour負(fù)二項(xiàng)分布\\'+str(vds)+' '+str(start)+'  '+' 符合負(fù)二項(xiàng)分布'+'  '+'.png')
    # if temp==2:
    #     plt.savefig('E:\\arrival\hour泊松分布\\'+str(vds)+' '+str(start)+'  '+' 符合泊松分布'+'  '+'.png')
    # if temp==3:
    #     plt.savefig('E:\\arrival\hour二項(xiàng)分布\\'+str(vds)+' '+str(start)+'  '+' 二項(xiàng)分布'+'  '+'.png')
    C.append(temp)
    A.append(vds)
    B.append(start)
    plt.close()




def huatu():#畫散點(diǎn)圖函數(shù)

    temp=PHT()
    plt.subplot(311)
    plt.scatter(K,V)
    plt.title(str(vds)+' K-V '+str(start))
    plt.subplot(312)
    plt.scatter(K,Q)
    plt.title(str(vds)+' K-Q '+str(start))
    plt.subplot(313)
    plt.title(str(vds)+' Q-V '+str(start))
    plt.scatter(Q,V)

    if temp == 0:
        plt.savefig('E:\\traffic\\'+str(vds)+' '+str(start)+' '+str(end)+'.png')
    else:
        plt.savefig('E:\\traffic\\' + str(vds) + '  ' + str(start) + ' ' + str(end) + str(temp) + '.png')
    plt.close()

vdss1=[1018,1026,1028,1030,1032,1035,1036,1016,1007]#向西方向
vdss2=[1019,1027,1029,1031,1033,1034,1037,1017,1008]#向東方向

# vdss1=[1008]
# for vds in vdss1:
#     for start in range(10,22):
#         init(start,start+1,vds)
#         hourarrival()
#


for vds in vdss1:
    for start in range(10,22):
        init(start,start+1,vds)
        hourarrival()
        #PHT()

for vds in vdss2:
    for start in range(10,22):
        init(start,start+1,vds)
        hourarrival()
        #PHT()


# arr=pd.DataFrame({"vds":A,'day':B,'M':M,'D':S,'kind':D,'time':C})
arr=pd.DataFrame({"vds":A,'day':B,'符合分布':C})
arr.to_excel('output1.xls')
# print(C)