MeteoInfoLab脚本示例：台风和热带气旋

qxsw2016

请教一个台风移动路径频数的问题，如有多年台风的移动路径数据，把区域分成等格距网格点，比如2.5度X2.5度，如何计算经过各网格点的次数，老师能不能给个思路。效果如图

MeteoInfo

qxsw2016 发表于 2018-9-13 21:41
请教一个台风移动路径频数的问题，如有多年台风的移动路径数据，把区域分成等格距网格点，比如2.5度X2.5度 ...

可以用topology模块的 intersects 函数判断两个几何图形是否相交，参考此脚本：

1. import mipylib.geolib.topology as tp
   
2. 
   
3. # Read typhoon data file
   
4. fn = 'D:/Temp/ascii/CH2015BST.txt'
   
5. tf = open(fn)
   
6. lons = []
   
7. lats = []
   
8. prss = []
   
9. for line in tf:
   
10.     #print line
    
11.     data = line.split()
    
12.     pn = int(data[2])
    
13.     for i in range(pn):
    
14.         line = tf.readline()
    
15.         data = line.split()
    
16.         lat = float(data[2])        
    
17.         lats.append(lat * 0.1)
    
18.         lon = float(data[3])
    
19.         lons.append(lon * 0.1)
    
20.         t = data[0]
    
21.         prs = float(data[4])
    
22.         prss.append(prs)
    
23.     lons.append(nan)
    
24.     lats.append(nan)
    
25.     prss.append(nan)
    
26. 
    
27. # Plot
    
28. axesm()
    
29. lworld = shaperead('D:/Temp/map/country1.shp')
    
30. geoshow(lworld, facecolor=[200,200,200])
    
31. layer = plotm(lons, lats, linewidth=1)
    
32. xlim(100, 210)
    
33. ylim(0, 60)
    
34. title('Typhoon pathway')
    
35. 
    
36. #Add a rectangle polygon
    
37. lat = array([20, 20, 10, 10, 20])
    
38. lon = array([180, 190, 190, 180, 180])
    
39. g1 = geoshow(lat, lon, displaytype='polygon', color=[150,230,230,230], edgecolor='b', size=2)
    
40. 
    
41. #Count the number of typhoons cross the rectangle
    
42. n = 0
    
43. for s in layer.shapes():
    
44.     if tp.intersects(s, g1):
    
45.         n += 1
    
46. print 'Number: %i' % n

import mipylib.geolib.topology as tp<br /> <br /> # Read typhoon data file<br /> fn = 'D:/Temp/ascii/CH2015BST.txt'<br /> tf = open(fn)<br /> lons = []<br /> lats = []<br /> prss = []<br /> for line in tf:<br />     #print line<br />     data = line.split()<br />     pn = int(data[2])<br />     for i in range(pn):<br />         line = tf.readline()<br />         data = line.split()<br />         lat = float(data[2])        <br />         lats.append(lat * 0.1)<br />         lon = float(data[3])<br />         lons.append(lon * 0.1)<br />         t = data[0]<br />         prs = float(data[4])<br />         prss.append(prs)<br />     lons.append(nan)<br />     lats.append(nan)<br />     prss.append(nan)<br /> <br /> # Plot<br /> axesm()<br /> lworld = shaperead('D:/Temp/map/country1.shp')<br /> geoshow(lworld, facecolor=[200,200,200])<br /> layer = plotm(lons, lats, linewidth=1)<br /> xlim(100, 210)<br /> ylim(0, 60)<br /> title('Typhoon pathway')<br /> <br /> #Add a rectangle polygon<br /> lat = array([20, 20, 10, 10, 20])<br /> lon = array([180, 190, 190, 180, 180])<br /> g1 = geoshow(lat, lon, displaytype='polygon', color=[150,230,230,230], edgecolor='b', size=2)<br /> <br /> #Count the number of typhoons cross the rectangle<br /> n = 0<br /> for s in layer.shapes():<br />     if tp.intersects(s, g1):<br />         n += 1<br /> print 'Number: %i' % n

qxsw2016

MeteoInfo 发表于 2018-9-13 23:18
可以用topology模块的 intersects 函数判断两个几何图形是否相交，参考此脚本：

好久没交流了，每次都有收获，感谢，我来学习学习！

qxsw2016

请教，如何把经过某个特定区域的台风路径挑选出来。即用topology模块的 intersects 函数判断两个几何图形是否相交，来判断是否经过某区域，但如何来标识区分这个台风呢？望百忙中指点。

MeteoInfo

qxsw2016 发表于 2018-9-17 22:18
请教，如何把经过某个特定区域的台风路径挑选出来。即用topology模块的 intersects 函数判断两个几何图形是 ...

可以给台风图层增加一个字段来记录是否经过某个区域。参考这里：http://bbs.06climate.com/forum.p ... &extra=page%3D1

qxsw2016

感谢老师 ，我来学习消化一下。还有一个问题，判断一个站点是否在某个行政区里，即点是否在某个区域里，是不也用 intersects 函数？

MeteoInfo

qxsw2016 发表于 2018-9-25 21:54
感谢老师 ，我来学习消化一下。还有一个问题，判断一个站点是否在某个行政区里，即点是否在某个区域里，是 ...

你可以试试，也可以参考这里：http://bbs.06climate.com/forum.p ... &extra=page%3D1

qxsw2016

学习了，太强大了，以前我做160个站是，硬是把每个城市所在省份找出来，