绘制气团轨迹动态图

MeteoInfo

HYSPLIT模式输出的气团轨迹文件可以通过MeteoInfoLab的addfile_hytraj函数读取，通过变量lon, lat, height读取气团轨迹经纬度和高度坐标。在地图坐标系中绘制地理底图，并用scatter函数绘制所有气团的终点位置（后向气团轨迹），然后通过循环绘制气团轨迹运动动画图。循环中的绘图函数plot来绘制气团轨迹线条，scatter来绘制该循环时次气团终点，需要注意的是下个循环绘制前用两个cll()函数删除上次循环绘制的轨迹线条和点。可以用time模块的sleep函数来控制动画的播放速度。

1. import time
   
2. 
   
3. fn = os.path.join(migl.get_sample_folder(), 'HYSPLIT', 'tdump')
   
4. f = addfile_hytraj(fn)
   
5. lon = f['lon'][:]
   
6. lat = f['lat'][:]
   
7. height = f['height'][:]
   
8. 
   
9. #Plot
   
10. ax = axesm()
    
11. geoshow('country', edgecolor='gray', facecolor=(230,230,230))
    
12. levs = arange(0, 600, 50)
    
13. scatter(lon[:,0], lat[:,0], size=6, color='r', marker='S')
    
14. yticks(arange(35, 54, 5))
    
15. title('Trajectory')
    
16. 
    
17. ntraj, npoint = lon.shape
    
18. for i in range(2, npoint+1):
    
19.     if i > 2:
    
20.         cll()
    
21.         cll()
    
22.     tlayer = plot(lon[:,-i:], lat[:,-i:], zvalues=height[:,-i:], levels=levs, linewidth=2)
    
23.     scatter(lon[:,-i], lat[:,-i], marker='s', size=5, edgecolor='m', facecolor=None)
    
24.     if i == 2:
    
25.         colorbar(tlayer)
    
26.     xlim(-92, -55)
    
27.     ylim(34, 54)
    
28.     plt.draw()
    
29.    
    
30.     time.sleep(0.1)

import time<br /> <br /> fn = os.path.join(migl.get_sample_folder(), 'HYSPLIT', 'tdump')<br /> f = addfile_hytraj(fn)<br /> lon = f['lon'][:]<br /> lat = f['lat'][:]<br /> height = f['height'][:]<br /> <br /> #Plot<br /> ax = axesm()<br /> geoshow('country', edgecolor='gray', facecolor=(230,230,230))<br /> levs = arange(0, 600, 50)<br /> scatter(lon[:,0], lat[:,0], size=6, color='r', marker='S')<br /> yticks(arange(35, 54, 5))<br /> title('Trajectory')<br /> <br /> ntraj, npoint = lon.shape<br /> for i in range(2, npoint+1):<br />     if i > 2:<br />         cll()<br />         cll()<br />     tlayer = plot(lon[:,-i:], lat[:,-i:], zvalues=height[:,-i:], levels=levs, linewidth=2)<br />     scatter(lon[:,-i], lat[:,-i], marker='s', size=5, edgecolor='m', facecolor=None)<br />     if i == 2:<br />         colorbar(tlayer)<br />     xlim(-92, -55)<br />     ylim(34, 54)<br />     plt.draw()<br />     <br />     time.sleep(0.1)

20231025_165640.mp4 (863.52 KB, 下载次数: 27)

2023-11-3 16:02 上传

点击文件名下载附件

将动画保存为gif动画图，用gifanimation函数生成一个动画对象，函数第一个参数是gif文件路径，第二个参数是gif动画每帧的切换时长。在循环中用gifaddframe函数将当前图像帧加入动画，最后用动画对象的finish函数结束动画gif图输出。

1. fn = os.path.join(migl.get_sample_folder(), 'HYSPLIT', 'tdump')
   
2. f = addfile_hytraj(fn)
   
3. lon = f['lon'][:]
   
4. lat = f['lat'][:]
   
5. height = f['height'][:]
   
6. 
   
7. #Plot
   
8. ax = axesm()
   
9. geoshow('country', edgecolor='gray', facecolor=(230,230,230))
   
10. levs = arange(0, 600, 50)
    
11. scatter(lon[:,0], lat[:,0], size=6, color='r', marker='S')
    
12. yticks(arange(35, 54, 5))
    
13. title('Trajectory')
    
14. 
    
15. giffn = 'D:/Temp/test/traj_animation.gif'
    
16. animation = gifanimation(giffn, delay=100)
    
17. ntraj, npoint = lon.shape
    
18. for i in range(2, npoint+1):
    
19.     print(i)
    
20.     if i > 2:
    
21.         cll()
    
22.         cll()
    
23.     tlayer = plot(lon[:,-i:], lat[:,-i:], zvalues=height[:,-i:], levels=levs, linewidth=2)
    
24.     scatter(lon[:,-i], lat[:,-i], marker='s', size=5, edgecolor='m', facecolor=None)
    
25.     if i == 2:
    
26.         colorbar(tlayer)
    
27.     xlim(-92, -55)
    
28.     ylim(34, 54)
    
29.     plt.draw()
    
30.     gifaddframe(animation)
    
31. 
    
32. animation.finish()

fn = os.path.join(migl.get_sample_folder(), 'HYSPLIT', 'tdump')<br /> f = addfile_hytraj(fn)<br /> lon = f['lon'][:]<br /> lat = f['lat'][:]<br /> height = f['height'][:]<br /> <br /> #Plot<br /> ax = axesm()<br /> geoshow('country', edgecolor='gray', facecolor=(230,230,230))<br /> levs = arange(0, 600, 50)<br /> scatter(lon[:,0], lat[:,0], size=6, color='r', marker='S')<br /> yticks(arange(35, 54, 5))<br /> title('Trajectory')<br /> <br /> giffn = 'D:/Temp/test/traj_animation.gif'<br /> animation = gifanimation(giffn, delay=100)<br /> ntraj, npoint = lon.shape<br /> for i in range(2, npoint+1):<br />     print(i)<br />     if i > 2:<br />         cll()<br />         cll()<br />     tlayer = plot(lon[:,-i:], lat[:,-i:], zvalues=height[:,-i:], levels=levs, linewidth=2)<br />     scatter(lon[:,-i], lat[:,-i], marker='s', size=5, edgecolor='m', facecolor=None)<br />     if i == 2:<br />         colorbar(tlayer)<br />     xlim(-92, -55)<br />     ylim(34, 54)<br />     plt.draw()<br />     gifaddframe(animation)<br /> <br /> animation.finish()

薛定谔的猫啊猫

请问下老师这个报错是什么意思啊indices must be 1 dimensions!

LLLLLYH

请问老师用plot3里画的时候怎么调轨迹的透明度呢？直接加alpha参数只能修改其中一个颜色的透明度

MeteoInfo

薛定谔的猫啊猫 发表于 2023-11-6 16:48
请问下老师这个报错是什么意思啊indices must be 1 dimensions!

可能是读取数组的维数设置问题

MeteoInfo

LLLLLYH 发表于 2023-11-6 17:32
请问老师用plot3里画的时候怎么调轨迹的透明度呢？直接加alpha参数只能修改其中一个颜色的透明度

这是二维图，没有用到plot3函数

LLLLLYH

对的，如果应用到plot3，请问有办法把轨迹透明度统一调整到半透明吗？

MeteoInfo

LLLLLYH 发表于 2023-11-7 14:59
对的，如果应用到plot3，请问有办法把轨迹透明度统一调整到半透明吗？

用alpha参数，比如  alpha=0.5

LLLLLYH

MeteoInfo 发表于 2023-11-7 18:19
用alpha参数，比如  alpha=0.5

我试过了，用在画HYSPLIT轨迹上的水汽含量时只能改变色标中的几个颜色
colsh2 = makecolors(len(levs) + 1, cmap='WhiteGreen', alpha=0.2)
traj = plot3(lon[ip,:], lat[ip,:], alt[ip,:], mvalues=sh, levels=levsh,colors=colsh2,linewidth=2)

MeteoInfo

LLLLLYH 发表于 2023-11-8 09:40
我试过了，用在画HYSPLIT轨迹上的水汽含量时只能改变色标中的几个颜色
colsh2 = makecolors(len(levs) + ...

应该是有bug，有空了我更新一下