Python处理云量信息

香草拿铁

期待作者的经验分享

junyang517

可以把代码贴上来看看，让大家给分析分析~~
另外，请问楼主的ECMWF数据是在哪里下载的？免费的吗？能否提供一下网址？

lwqqqqq

chongzika 发表于 2018-12-9 10:38
回头楼主可以分享下经验

尴尬...就用的for循环，处理一个月的数据大概需要20分钟 很慢

lwqqqqq

香草拿铁 发表于 2018-12-9 20:47
期待作者的经验分享

没有经验分享 用的for循环，尝试了一次多线程，还不如单线程处理呢

lwqqqqq

junyang517 发表于 2018-12-10 12:39
可以把代码贴上来看看，让大家给分析分析~~
另外，请问楼主的ECMWF数据是在哪里下载的？免费的吗？能否提 ...

对的 是免费的 网址如下：https://www.ecmwf.int/forecasts/datasets/archive-datasets/browse-reanalysis-datasets
代码稍后贴上

lwqqqqq

import threading
from netCDF4 import Dataset
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap
import time


def data_process(month):
    nc_obj = Dataset("E:\\data\\2015_" + month + ".nc")
   
    month_dict={'Jan':1,'Feb':2,'Mar':3,'Apr':4,'May':5,'Jun':6,
                'Jul':7,'Aug':8,'Sep':9,'Oct':10,'Nov':11,'Dec':12}
   
    v = [month_dict[k] for k in month_dict.keys() if k==month][0]
    list1 = [1,3,5,7,8,10,12]
    list2 = [2]
    list3 = [4,6,9,11]
        
    lats_list = range(53)
   
    lons_list = range(81)
   
    if v in list1:
        times_list = range(8*31)
    if v in list2:
        times_list = range(8*28)
    if v in list3:
        times_list = range(8*30)
        
    arr1 = np.zeros( (53,81), dtype=np.float64 )

    for lat in lats_list:
            for lon in lons_list:
                for t in times_list:
                    hcc = nc_obj.variables["hcc"][t,lat,lon]
                    lcc = nc_obj.variables["lcc"][t,lat,lon]
                    mcc = nc_obj.variables["mcc"][t,lat,lon]
                    count = 0
                    if lcc <= 0.05 and mcc <= 0.05 and hcc >= 0.05:
                        count += 1
                arr1[lat][lon] = count/len(times_list)


    m=Basemap(projection='tmerc',lat_0 = 23, lon_0 = 118,llcrnrlat=16.5,
              llcrnrlon=108,urcrnrlat=23,urcrnrlon=118)
#加载中国地图
    m.readshapefile('E:\\China_map\\region','region.shp',color='red')
   
#生成绘图背景格点场
    ny=arr1.shape[0];nx=arr1.shape[1]
    lons,lats=m.makegrid(nx,ny)
    x,y=m(lons,lats)
    yy=y[::-1]
   
#画经纬度网格
    m.drawparallels(np.arange(23,16.5,-1.5),labels=[1,0,0,0])
    m.drawmeridians(np.arange(108,118,1.5),labels=[0,0,0,1])
   
#画填色图
    cm = plt.cm.get_cmap('Blues')    #读取一个内置的色系
    shades=m.contourf(x,yy,arr1,cmap=cm)  #将数据填色绘出
    m.colorbar(shades)                      #添加colorbar

#画等值线
#curve=m.contour(x,yy,arr1,colors='black')
#plt.clabel(curve,fmt='%1.0f',colors='black')
   
#添加标题
    plt.title('2015_%d_hcc'%v)
   
#保存
    plt.savefig('2015_hcc\\%d月.png'%v,dpi=300)
#关闭绘图功能
    plt.close()


   
#单线程处理
#time1 = time.time()
#data_process("Jun")
#time2 = time.time()
#print(time2 - time1)
#print("\a")


#多线程处理
threads = []
t1 = threading.Thread(target=data_process,args=('Sep',))
threads.append(t1)
t2 = threading.Thread(target=data_process,args=('Nov',))
threads.append(t2)

time1 = time.time()

if __name__ == '__main__':
    for t in threads:
        t.setDaemon(True)
        t.start()
    t.join()
    print("\a")

time2 = time.time()
print(time2 - time1)

lwqqqqq

lwqqqqq 发表于 2018-12-10 16:34
对的 是免费的 网址如下：https://www.ecmwf.int/forecasts/datasets/archive-datasets/browse-reanalysi ...

代码已发，麻烦给点意见，看看能不能优化下

暗藏

tmp = np.zeros((len(times_list),53,81), dtype=np.float64)
hcc = nc_obj.variables["hcc"]   # 不清楚你的文件结构，这部分酌情替换为nc数据读入一个三维数组
lcc = nc_obj.variables["lcc"]
mcc = nc_obj.variables["mcc"]
tmp[(lcc <= 0.05)&(mcc <= 0.05)&(hcc >= 0.05)] = 1
arr1 = np.mean(tmp,0)

以上替换for的部分，我没有测试，仅供参考

lwqqqqq

暗藏 发表于 2018-12-12 11:28
tmp = np.zeros((len(times_list),53,81), dtype=np.float64)
hcc = nc_obj.variables["hcc"]   # 不清楚 ...

好的 我来试一下 多谢了

lwqqqqq

暗藏 发表于 2018-12-12 11:28
tmp = np.zeros((len(times_list),53,81), dtype=np.float64)
hcc = nc_obj.variables["hcc"]   # 不清楚 ...

太厉害了！！！用for循环得20分钟，改了一下只要3秒，为什么这么快呢？太厉害了