- 积分
- 10278
- 贡献
-
- 精华
- 在线时间
- 小时
- 注册时间
- 2019-7-29
- 最后登录
- 1970-1-1
|
登录后查看更多精彩内容~
您需要 登录 才可以下载或查看,没有帐号?立即注册
x
本帖最后由 雨落森林 于 2023-4-17 18:28 编辑
代码改自姚素香老师的Fortran程序,不过由于咱们使用Python,不可能用好几层for循环来套,咱们用矩阵思想,对格点求偏导,如T[...,1:-1,:],即可指定空间的纬度方向上第2个到倒数第2个格点的温度,以此类推。思路和公式见docx文档,代码我也写好了注释。数据来自NCEP,是nc格式,经过提取为npy,其中,T、U、V等变量的shape是41,282,17,73,144,即41年,248个时次(7和8月、逐6小时),17个高度层,73个纬度(-90到90,间隔2.5),144个经度(0到357.5,间隔2.5),由于数据太大占用硬盘空间28G,计算时需要四五十个G的内存,故本帖只放代码,不放数据了,仅供学习一个思路。最后见一下成果图吧,图给的是某年7、8月平均的1000hPa的非绝热加热,单位是(K/s),可以看到,在青藏高原有个特别大的加热,说明算出来是对的。附件里面的压缩文件有三个代码,data_process.py是数据的提取和处理(反转纬度),主要看main.py了解具体思路(每项都计算后存储在part变量里面,每次Q加上part后清零),test.py是用来画图。最后需要注意的是,计算出来的量级应当是10的-5次方的量级,否则就有可能是哪里错了。 |
|