python 您所在的位置:网站首页 怎样读地形图的经纬度 python

python

2024-07-04 15:33| 来源: 网络整理| 查看: 265

0.写在前面

本文接着上一篇博客,将用imshow函数把数目分布叠加在背景上,不再需要用其他软件把图片叠起来了。

1.话不多说上代码

(上篇文章中部分篇幅写得更详细些,包括文件读取等,也可作为参考)

import matplotlib.pyplot as plt import pylab import numpy as np import xarray as xr from cartopy.mpl.ticker import LongitudeFormatter, LatitudeFormatter import cartopy.crs as ccrs import pandas as pd from cartopy.io import shapereader as shpreader from cartopy.io.shapereader import Reader plt.rcParams['font.sans-serif'] = ['SimHei'] plt.rcParams['axes.unicode_minus'] = False df = pd.read_table(r'D:\何卓尧\桌面\2014-2020\2014\SURF_CHN_MUL_DAY_20140101.txt',sep='\t',header=None) #下面是站点经纬度数据 lat = df[64] del lat[0] lat=lat.astype(float) lon = df[66] del lon[0] lon=lon.astype(float) ds = xr.open_dataset(r"地形图\ETOPO2v2c_f4.nc") # 下面是准备用于绘制地形图的数据 lon1 = np.linspace(min(ds['x'].data), max( ds['x'].data), len(ds['x'].data)) # 经度 lat1 = np.linspace(min(ds['y'].data), max( ds['y'].data), len(ds['y'].data)) # 纬度 dem = ds['z'].data #准备画布 fig = plt.figure(1, figsize=[16, 9]) proj = ccrs.PlateCarree() ax = plt.subplot(1, 1, 1, projection=proj) extent = [70, 140, 30, 55] ax.set(xlim=(70,140),ylim=(30,55)) #用numb数组保存每个格子的站点数量 numb=np.zeros([50,140],dtype=float,order='C') for m in range(1,1740): if lat[m]=30 and lon[m]>=70 and lon[m]


【本文地址】

公司简介

联系我们

今日新闻

    推荐新闻

    专题文章
      CopyRight 2018-2019 实验室设备网 版权所有