如何在柱坐标系中进行插值 - 编程论坛

RT，小弟想要实现的功能：已知一个圆柱形建筑物的几个沉降观测点（坐标X、Y）的沉降数据（Z），绘制一个基础的沉降量云图
一开始直接用X、Y插值，但是绘制出来的图是多边形（我用的是12点的数据，要是只有4个点数据，就变成四边行了），而不是圆形的。 →  interpolateb_in_cat
后来改用柱坐标系，先生成r、theta的meshgrid，再转换成x、y插值，但是还是多边形，始终无法绘制出一个完整的圆形来。  →  interpolate_in_pol
又尝试了柱坐标系的r、theta直接插值，再转成x、y，但是啥结果也没有了……  →  interpolate_in_pol2
请教各位大神，如何在柱坐标系中进行插值，并绘制成一个完整的圆形云图？谢了先

并祝新年快乐！

只有本站会员才能查看附件，请登录

程序代码：

import matplotlib.pyplot as plt
import numpy as np
from scipy.interpolate import griddata
R=10
# 12 Point Sample Data
x = np.array([5,8.660254038,10,8.660254038,5,1.22515E-15,-5,-8.660254038,-10,-8.660254038,-5,-2.4503E-15])
y = np.array([8.660254038,5,6.12574E-16,-5,-8.660254038,-10,-8.660254038,-5,-1.83772E-15,5,8.660254038,10])
z = np.array([0.006,0.021,0.063,0.028,0.012,0.015,0.119,0.07,0.071,0,0.003,0.054])

def pol2cart(r, theta):
    x = r* np.cos(theta)
    y = r* np.sin(theta)
    return (x, y)

def cart2pol(x, y):
    r = (x**2+y**2)**0.5
    theta = np.arctan2(y,x)
    return (r, theta)

# X = np.linspace(min(x), max(x))
# Y = np.linspace(min(y), max(y))
# X=r*np.cos(theta)
# Y=r*np.sin(theta)
# X, Y = np.meshgrid(X, Y)

def interpolate_in_pol(x,y,z,resolution = 50,contour_method='linear'):
    r = np.linspace(0, R,resolution)
    theta = np.linspace(0, 2 * np.pi,resolution)
    r, theta = np.meshgrid(r, theta)
    X, Y = pol2cart(r, theta)
    points = [[a,b] for a,b in zip(x,y)]
    Z = griddata(points, z, (X, Y), method=contour_method)
    return X,Y,Z

def interpolate_in_pol2(x,y,z,resolution = 50,contour_method='linear'):
    r,theta=cart2pol(x, y)
    R = np.linspace(0, r,resolution)
    Theta = np.linspace(0, 2 * np.pi,resolution)
    R, Theta = np.meshgrid(R, Theta)
    points = [[a,b] for a,b in zip(r,theta)]
    Z = griddata(points, z, ( R, Theta), method=contour_method)
    X, Y = pol2cart(R, Theta)
    return X,Y,Z

def interpolateb_in_cat(x,y,z,resolution = 50,contour_method='linear'):
    resolution = str(resolution)+'j'
    X,Y = np.mgrid[min(x):max(x):complex(resolution),   min(y):max(y):complex(resolution)]
    points = [[a,b] for a,b in zip(x,y)]
    Z = griddata(points, z, (X, Y), method=contour_method)
    return X,Y,Z

X1,Y1,Z1 = interpolateb_in_cat(x,y,z,resolution = 500,contour_method='cubic')
X2,Y2,Z2  =  interpolate_in_pol(x,y,z,resolution = 500,contour_method='cubic')
X3,Y3,Z3 = interpolate_in_pol(x,y,z,resolution = 500,contour_method='cubic')

plt.figure()
plt.pcolormesh(X1, Y1, Z1)
plt.colorbar() # Color Bar
plt.show()

plt.figure()
plt.pcolormesh(X2, Y2, Z2)
plt.colorbar() # Color Bar
plt.show()

[此贴子已经被作者于2021-1-3 00:15编辑过]