Python实现绘制Matlab格式的地图边框的示例代码

2024-04-02 19:55

短信预约 -IT技能 免费直播动态提醒

1、Python绘制色斑图


import matplotlib.pyplot as plt
import numpy as np
from cartopy.mpl.ticker import LongitudeFormatter, LatitudeFormatter
import cartopy.crs as ccrs
import cartopy.feature as cfeature
import pandas as pd
import maskout
import matplotlib.patches as mpatches
from cartopy.io.shapereader import Reader
from matplotlib import rcParams
config = {"font.family":'Times New Roman',"font.size": 16,"mathtext.fontset":'stix'}
rcParams.update(config)
df1 = pd.read_excel(r"F:/Rpython/lp37/henanmap/henan.xlsx")
lat = df1['lat']
lon = df1['lon']
PM25 = df1['PM25']
# 设置经纬度并用scipy包中的函数插值
from scipy.interpolate import Rbf
olon = np.linspace(110, 117, 100)
olat = np.linspace(31, 37, 100)
olon, olat = np.meshgrid(olon, olat)
# 插值处理,‘linear'，‘nearest'，‘cubic'
# cubic, gaussian, inverse_multiquadric, linear, multiquadric, quintic, thin_plate
# rain_data_new = griddata((lon,lat), data, (olon,olat), method='linear')
func1=Rbf(lon,lat,PM25,function='linear')
PM25=func1(olon,olat)
print(olon)
print(olat)
print(PM25)
print(olat.shape)
print(PM25.shape)
# 建立画布
fig2 = plt.figure(figsize=(16, 12))
proj = ccrs.PlateCarree()
leftlon, rightlon, lowerlat, upperlat = (110, 117, 31, 37)  # 根据上下限确定范围，至少为10°
# 春 ------------------------------------------------------------------------------------
#左 底 宽 高
ax=fig2.add_subplot(1,1,1,projection=ccrs.PlateCarree())
# 在画布的绝对坐标建立子图
ax.set_extent([leftlon, rightlon, lowerlat, upperlat],crs=ccrs.PlateCarree())
# 设置地图属性:加载国界、海岸线，land 为灰色
ax.stock_img()
# 绘制河南省行政边界
ticks=np.arange(30,80,5)
cf=ax.contourf(olon,olat,PM25,levels=ticks,cmap='gist_rainbow',transform=ccrs.PlateCarree(),extend='both')
clip1=maskout.shp2clip(cf,ax,'F:/Rpython/lp37/henanmap/henan1')
cf=plt.colorbar(cf,ticks=ticks,shrink=0.96,orientation='vertical',extend='both',pad=0.01,aspect=35)
ax.add_feature(cfeature.COASTLINE.with_scale('50m'),linewidth=0.5,zorder=2,color='k')# 添加海岸线
ax.add_feature(cfeature.LAKES.with_scale('50m'))
ax.add_feature(cfeature.RIVERS.with_scale('50m'))
ax.add_feature(cfeature.OCEAN.with_scale('50m'))
ax.add_feature(cfeature.LAND.with_scale('50m'))
ax.add_geometries(Reader(r'F:/Rpython/lp37/henanmap/henan1.shp').geometries(),ccrs.PlateCarree(),facecolor='none',edgecolor='k',linewidth=1.2,zorder=1)
ax.add_geometries(Reader(r'F:/Rpython/lp27/data/river1.shp').geometries(),ccrs.PlateCarree(),facecolor='none',edgecolor='b',linewidth=0.2)
ax.add_geometries(Reader(r'F:/Rpython/lp27/data/china1.shp').geometries(),ccrs.PlateCarree(),facecolor='none',edgecolor='r',linewidth=1.2)
ax.add_geometries(Reader(r'F:/Rpython/lp27/data/china2.shp').geometries(),ccrs.PlateCarree(),facecolor='none',edgecolor='r',linewidth=0.8)
ax.add_geometries(Reader(r'F:/Rpython/lp27/data/ne_50m_lakes.shp').geometries(),ccrs.PlateCarree(),facecolor='none',edgecolor='k',linewidth=0.2)
ax.add_geometries(Reader(r'F:/Rpython/lp37/data37/river/1级河流.shp').geometries(),ccrs.PlateCarree(),facecolor='none',edgecolor='RoyalBlue',linewidth=0.4)
ax.add_geometries(Reader(r'F:/Rpython/lp37/data37/river/2级河流.shp').geometries(),ccrs.PlateCarree(),facecolor='none',edgecolor='DodgerBlue',linewidth=0.3)
ax.add_geometries(Reader(r'F:/Rpython/lp37/data37/river/3级河流.shp').geometries(),ccrs.PlateCarree(),facecolor='none',edgecolor='DeepSkyBlue',linewidth=0.2)
ax.add_geometries(Reader(r'F:/Rpython/lp37/data37/river/4级河流.shp').geometries(),ccrs.PlateCarree(),facecolor='none',edgecolor='SkyBlue',linewidth=0.15)
ax.add_geometries(Reader(r'F:/Rpython/lp37/data37/river/5级河流.shp').geometries(),ccrs.PlateCarree(),facecolor='none',edgecolor='LightSkyBlue',linewidth=0.05)
ax.add_geometries(Reader(r'F:/Rpython/lp37/data37/river/主要湖泊.shp').geometries(),ccrs.PlateCarree(),edgecolor='none',linewidth=0,facecolor='#BEE8FF')
# 以下6条语句是定义地理坐标标签格式
ax.set_xticks(np.arange(leftlon, rightlon+0.1, 1),crs=ccrs.PlateCarree())
ax.set_yticks(np.arange(lowerlat, upperlat+0.1, 1),crs=ccrs.PlateCarree())
lon_formatter = LongitudeFormatter()
lat_formatter = LatitudeFormatter()
ax.xaxis.set_major_formatter(lon_formatter)
ax.yaxis.set_major_formatter(lat_formatter)
plt.savefig('F:/Rpython/lp37/plot185.4.png',dpi=600,bbox_inches='tight',pad_inches=0)
plt.show()

2、Python绘制比例尺、指南针


# 添加比例尺，指南针
def add_north(ax, labelsize=18, loc_x=0.95, loc_y=0.99, width=0.06, height=0.09, pad=0.14):
    """
    画一个比例尺带'N'文字注释
    主要参数如下
    :param ax: 要画的坐标区域 Axes实例 plt.gca()获取即可
    :param labelsize: 显示'N'文字的大小
    :param loc_x: 以文字下部为中心的占整个ax横向比例
    :param loc_y: 以文字下部为中心的占整个ax纵向比例
    :param width: 指南针占ax比例宽度
    :param height: 指南针占ax比例高度
    :param pad: 文字符号占ax比例间隙
    :return: None
    """
    minx, maxx = ax.get_xlim()
    miny, maxy = ax.get_ylim()
    ylen = maxy - miny
    xlen = maxx - minx
    left = [minx + xlen*(loc_x - width*.5), miny + ylen*(loc_y - pad)]
    right = [minx + xlen*(loc_x + width*.5), miny + ylen*(loc_y - pad)]
    top = [minx + xlen*loc_x, miny + ylen*(loc_y - pad + height)]
    center = [minx + xlen*loc_x, left[1] + (top[1] - left[1])*.4]
    triangle = mpatches.Polygon([left, top, right, center], color='k')
    ax.text(s='N',
            x=minx + xlen*loc_x,
            y=miny + ylen*(loc_y - pad + height),
            fontsize=labelsize,
            horizontalalignment='center',
            verticalalignment='bottom')
    ax.add_patch(triangle)
#-----------函数：添加比例尺--------------
def add_scalebar(ax,lon0,lat0,length,size=0.01):
    '''
    ax: 坐标轴
    lon0: 经度
    lat0: 纬度
    length: 长度
    size: 控制粗细和距离的
    '''
    # style 3
    ax.hlines(y=lat0,  xmin = lon0, xmax = lon0+length/111, colors="black", ls="-", lw=1,linewidth=3, label='%d km' % (length))
    ax.vlines(x = lon0, ymin = lat0-size, ymax = lat0+size, colors="black", ls="-", lw=1,linewidth=3)
    ax.vlines(x = lon0+length/2/111, ymin = lat0-size, ymax = lat0+size, colors="black", ls="-", lw=1,linewidth=3)
    ax.vlines(x = lon0+length/111, ymin = lat0-size, ymax = lat0+size, colors="black", ls="-", lw=1,linewidth=3)
    ax.text(lon0+length/111,lat0+size+0.05,'%d' % (length),horizontalalignment = 'center')
    ax.text(lon0+length/2/111,lat0+size+0.05,'%d' % (length/2),horizontalalignment = 'center')
    ax.text(lon0,lat0+size+0.05,'0',horizontalalignment = 'center')
    ax.text(112.28,31.45,'km',horizontalalignment = 'center')

3、Python绘制Matlab格式的地图边框


def drow_the_scale(y,x,text,length = 1.5,lw = 5):
    #画比例尺函数
    # y代表比例尺所在纬度
    # x代表比例尺开始的经度
    # text代表比例尺最后刻度值
    # length代表比例尺的长度，单位为多少个经度
    # lw代表比例尺的宽度
    step = length/5#计算步长，画五格
    #画黑白线五条
    plt.hlines(y=y,xmin=x,xmax=x + step,colors="black", ls="-", lw=lw)
    plt.hlines(y=y,xmin=x + step,xmax=x + step*2,colors="white", ls="-", lw=lw)
    plt.hlines(y=y,xmin=x + step*2,xmax=x + step*3,colors="black", ls="-", lw=lw)
    plt.hlines(y=y,xmin=x + step*3,xmax=x + step*4,colors="white", ls="-", lw=lw)
    plt.hlines(y=y,xmin=x + step*4,xmax=x + step*5,colors="black", ls="-", lw=lw)
    #画长刻度两个
    plt.vlines(x = x, ymin = y - (lw/100) *3, ymax = y + lw/100, colors="black", ls="-", lw=1)
    plt.vlines(x = x + length, ymin = y - (lw/100) *3, ymax = y + lw/100, colors="black", ls="-", lw=1)
    #画段刻度四个
    plt.vlines(x = x + step, ymin = y - (lw/100) *2, ymax = y + lw/100, colors="black", ls="-", lw=1)
    plt.vlines(x = x + step*2, ymin = y - (lw/100) *2, ymax = y + lw/100, colors="black", ls="-", lw=1)
    plt.vlines(x = x + step*3, ymin = y - (lw/100) *2, ymax = y + lw/100, colors="black", ls="-", lw=1)
    plt.vlines(x = x + step*4, ymin = y - (lw/100) *2, ymax = y + lw/100, colors="black", ls="-", lw=1)
    #写字，0，500，km
    plt.text(x,y - (lw/100) *7,'0',horizontalalignment = 'center')
    plt.text(x + length,y - (lw/100) *7,text,horizontalalignment = 'center')
    plt.text(x + length/2,y + (lw/100)*2,'km',horizontalalignment = 'center')
def drowscale(extent,scale_y,scale_x,scale_text,step = 5,lw = 10,scale_length = 1.5,scale_lw = 5):
    # 画地图黑白边框和比例尺
    # extent：表示四周经纬度[west, east, south, north]
    # scale_y,scale_x,scale_text：代表比例尺的位置，纬度，经度，刻度值
    # step：表示步长，一格代表几个经纬度
    # lw：代表边框宽度
    # scale_length：代表比例尺长度（单位为经度例如1.5个经度）
    # scale_lw：代表比例尺宽度
    for y in [extent[2],extent[3]] :#画上下两边框
        xmin = extent[0]
        while (xmin < extent[1]):
            plt.hlines(y=y,xmin=xmin,xmax=xmin+step,colors="white", ls="-", lw=lw)
            xmin = xmin+step*2
        xmin = extent[0]+step
        while (xmin < extent[1]):
            plt.hlines(y=y,xmin=xmin,xmax=xmin+step,colors="black", ls="-", lw=lw)
            xmin = xmin+step*2
    for x in [extent[0],extent[1]] :#画左右两边狂
        ymin = extent[2]
        while (ymin < extent[3]):
            plt.vlines(x = x, ymin = ymin, ymax = ymin+step, colors="black", ls="-", lw=lw)
            ymin = ymin+step*2
        ymin = extent[2]+step
        while (ymin < extent[3]):
            plt.vlines(x = x, ymin = ymin, ymax = ymin+step, colors="white", ls="-", lw=lw)
            ymin = ymin+step*2
    drow_the_scale(scale_y,scale_x,scale_text)#画比例尺