Python完美白化

猫老师007

大佬，有基于cartopy的白化教程吗

言墨

猫老师007 发表于 2020-4-29 15:53
大佬，有基于cartopy的白化教程吗

同问  我也用的cartopy  等经纬度投影可以用  但是其余的好像不太可以

洗雨

import shapefile
from matplotlib.path import Path
from matplotlib.patches import PathPatch
import requests
import json
import os
import datetime as DT
import time
from datetime import datetime
# 引用部分
import numpy as np
import pandas as pd
from scipy.interpolate import Rbf  # 径向基函数 ： 将站点信息插到格点上 用于绘制等值线
import test
import matplotlib.pyplot as plt
import matplotlib.colors as colors
import matplotlib as mpl
from mpl_toolkits.basemap import Basemap
import cmaps


dt1=DT.datetime(2015,5,2,0,0,0)
dt2=DT.datetime(2015,5,2,0,0,0)
utcDt1 = dt1
utcDt2 = dt2
cur0 = pd.date_range(utcDt1,utcDt2,freq='H')
def basemask(cs, ax, map, shpfile):

    sf = shapefile.Reader(shpfile)
    vertices = []
    codes = []
    for shape_rec in sf.shapeRecords():
        # if shape_rec.record[4] in shpfile:
        if shape_rec.record[0] >= 0:  
            pts = shape_rec.shape.points
            prt = list(shape_rec.shape.parts) + [len(pts)]
            for i in range(len(prt) - 1):
                for j in range(prt, prt[i+1]):
                    vertices.append(map(pts[j][0], pts[j][1]))
                codes += [Path.MOVETO]
                codes += [Path.LINETO] * (prt[i+1] - prt -2)
                codes += [Path.CLOSEPOLY]
            clip = Path(vertices, codes)
            clip = PathPatch(clip, transform = ax.transData)   
    for contour in cs.collections:
        contour.set_clip_path(clip)
        print('1234')
    return clip

for time1 in cur0:
    fn = time1.strftime('%Y%m%d%H')
    fid = r'E:\precipitation_test\hour_1\{0}.txt'.format(fn)
    #print(fid = input("请输入需要画图的时次和文本格式（r'E:\precipitation_test\hour_1\2015050100.txt')："))
    name = ['station','lon','lat','PRE_1h','TEM']
    data = pd.read_csv(fid,sep='\s+',header=None,names=name)
    # print(data)

    ######插值#########
    lon = data['lon']
    lat = data['lat']
    rain_data = data['PRE_1h']
    tem_data  = data['TEM']
    olon = np.linspace(78,100,88)
    olat = np.linspace(26,38,88)
    olon,olat = np.meshgrid(olon,olat)

    ######插值处理######
    func = Rbf(lon,lat,tem_data,function = 'linear')
    tem_data_new = func(olon,olat)

    plt.rcParams['font.sans-serif']=['SimHei'] #用来正常显示中文
    plt.rcParams['axes.unicode_minus']=False #用来正常显示负号

    #画布及绘图声明
    fig = plt.figure(figsize = (16,9.6))  # 画布
    # plt.rc('font',size=15,weight='bold')
    ax = fig.add_subplot(111) #绘图区
    m = Basemap(projection='cyl',llcrnrlat=26,llcrnrlon=78,urcrnrlat=38,urcrnrlon=100)
    m.readshapefile('E:/precipitation_test/Tibet_shps/tibet_shp/xizang_all','xizang_all.shp',default_encoding='gbk',drawbounds=True,linewidth=1, color='k')
    m.readshapefile('E:/precipitation_test/Tibet_shps/tibet_shp/river_1','river_1.shp',default_encoding='gbk',linewidth=1,color ='b')
    m.readshapefile('E:/precipitation_test/Tibet_shps/tibet_shp/river_2','river_2.shp',default_encoding='gbk',linewidth=0.8,color ='b')
    m.readshapefile('E:/precipitation_test/Tibet_shps/tibet_shp/river_3','river_3.shp',default_encoding='gbk',linewidth=0.6,color ='b')
    m.readshapefile('E:/precipitation_test/Tibet_shps/tibet_shp/xzlake','xzlake.shp',default_encoding='gbk',linewidth=1,color ='b')
    x,y = m(olon,olat)
    xx,yy = m(lon,lat)
    levels = np.linspace(0,np.max(tem_data_new),50)



    #色彩定制：24小时降水量级色标
    clevs = [-15.,-10.,-7.,-4.,-2.0,0.,2.,4.,7.,10.,15.] #自定义颜色列表


    #  绘制等值线、等值线填色
    ct = ax.contourf(x,y,tem_data_new,levels=clevs,cmap='seismic')
    plt.title('Tibetan Distribution of Temperature')
    position = fig.add_axes([0.91,0.15,0.01,0.7]) #位置【左，下，宽,高】


    plt.colorbar(ct,cax=position,label='Temperature (℃)')           #颜色参照表

    # 白化

    basemask(ct, ax, m, 'southwest')  
    print('hahaha')      




runfile('E:/precipitation_test/tem_test.py', wdir='E:/precipitation_test')
Traceback (most recent call last):


  File "E:\precipitation_test\tem_test.py", line 114, in <module>
    basemask(ct, ax, m, 'southwest')


  File "E:\precipitation_test\tem_test.py", line 32, in basemask
    for shape_rec in sf.shapeRecords():


  File "D:\ProgramData\Annaconda3\lib\site-packages\shapefile.py", line 1039, in shapeRecords
    for rec in zip(self.shapes(), self.records())])


  File "D:\ProgramData\Annaconda3\lib\site-packages\shapefile.py", line 1012, in records
    r = self.__record(oid=i)


  File "D:\ProgramData\Annaconda3\lib\site-packages\shapefile.py", line 987, in __record
    value = u(value, self.encoding, self.encodingErrors)


  File "D:\ProgramData\Annaconda3\lib\site-packages\shapefile.py", line 104, in u
    return v.decode(encoding, encodingErrors)


UnicodeDecodeError: 'utf-8' codec can't decode byte 0xc7 in position 0: invalid continuation byte


raceback (most recent call last):

  File "E:\precipitation_test\tem_test.py", line 18, in <module>
    import maskout

ValueError: source code string cannot contain null bytes


runfile('E:/precipitation_test/tem_test.py', wdir='E:/precipitation_test')
Traceback (most recent call last):

  File "E:\precipitation_test\tem_test.py", line 18, in <module>
    import maskout

ModuleNotFoundError: No module named 'maskout'
  运行clip = maskout.shp2clip(ct,ax,'E:\\precipitation_test\\Tibet_shps\\bou2_4p',[540000])报错应该怎么解决，我把utf-8改了还是报错

花生阿盛

学习，谢谢分享

洗雨

lovechang1314 发表于 2016-3-22 17:02
简单测试了以下，可以用。
有个问题，projection='cyl' 应该是不需要转换坐标的吧？
如果用'stere'这个时 ...

你们用的Python版本是哪个版本？我用的3.7报错
odeDecodeError: 'utf-8' codec can't decode byte 0xce in position 0: invalid continuation byte

洗雨

您们用的Python版本是哪个版本？我用的3.7报错
return v.decode(encoding, encodingErrors)
UnicodeDecodeError: 'utf-8' codec can't decode byte 0xce in position 0: invalid continuation byte

修改maskout 文件下utf-8为gbk或者gb2312或者utf-16都不行

平流层的萝卜

洗雨 发表于 2020-5-18 17:41
您们用的Python版本是哪个版本？我用的3.7报错
return v.decode(encoding, encodingErrors)
UnicodeDecod ...

这个是2.7的版本。你错误的原因有没有更详细的说明？

洗雨

平流层的萝卜 发表于 2020-5-18 19:55
这个是2.7的版本。你错误的原因有没有更详细的说明？

runfile('E:/precipitation_test/tem_test.py', wdir='E:/precipitation_test')
Reloaded modules: maskout01
Traceback (most recent call last):

  File "E:\precipitation_test\tem_test.py", line 116, in <module>
    basemask(ct, ax, m, 'Tibet')

  File "E:\precipitation_test\tem_test.py", line 33, in basemask
    for shape_rec in sf.shapeRecords():

  File "D:\ProgramData\Annaconda3\lib\site-packages\shapefile.py", line 1039, in shapeRecords
    for rec in zip(self.shapes(), self.records())])

  File "D:\ProgramData\Annaconda3\lib\site-packages\shapefile.py", line 1012, in records
    r = self.__record(oid=i)

  File "D:\ProgramData\Annaconda3\lib\site-packages\shapefile.py", line 987, in __record
    value = u(value, self.encoding, self.encodingErrors)

  File "D:\ProgramData\Annaconda3\lib\site-packages\shapefile.py", line 104, in u
    return v.decode(encoding, encodingErrors)

UnicodeDecodeError: 'utf-8' codec can't decode byte 0xce in position 0: invalid continuation byte

平流层的萝卜

洗雨 发表于 2020-5-19 12:06
runfile('E:/precipitation_test/tem_test.py', wdir='E:/precipitation_test')
Reloaded modules: mas ...

maskout是我编的，当时只考虑2.7的，没考虑3以上。但你提到3.5可以成功，结合你的错误信息，我觉得可能是你版本的shapefile模块有问题，导致编码错误，建议update一下pyshp模块和其他的相关的处理shapefile的模块。

洗雨

平流层的萝卜 发表于 2020-5-19 20:07
maskout是我编的，当时只考虑2.7的，没考虑3以上。但你提到3.5可以成功，结合你的错误信息，我觉得可能是 ...

好的，谢谢