【问题标题】:How to draw random planes如何绘制随机平面
【发布时间】:2014-10-02 20:33:34
【问题描述】:

我正在使用以下代码在 3d 中绘制穿过原点的随机平面。

from __future__ import division
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

#Number of hyperplanes
n = 20
#Dimension of space
d = 3

plt3d = plt.figure().gca(projection='3d')
for i in xrange(n):
    #Create random point on unit sphere
    v = np.random.normal(size = d)
    v = v/np.sqrt(np.sum(v**2))
    # create x,y
    xx, yy = np.meshgrid(range(-5,5), range(-5,5))
    z = (-v[0] * xx - v[1] * yy)/v[2]
    # plot the surface
    plt3d.plot_surface(xx, yy, z, alpha = 0.5)
plt.show()

但是看图片我不相信他们是统一选择的。我做错了什么?

【问题讨论】:

标签: python math numpy


【解决方案1】:

您的代码正在生成具有随机分布法线的平面。它们只是看起来不那样,因为 z 尺度比 x 和 y 尺度大得多。

您可以通过生成点来生成更好看的图像 均匀分布在平面上。为此,将平面参数化为 新坐标(u,v),然后在均匀间隔的网格上对平面进行采样 (u,v) 点。然后将这些 (u,v) 点转换为 (x,y,z)-空间中的点。

from __future__ import division
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
import math
import itertools as IT

def points_on_sphere(dim, N, norm=np.random.normal):
    """
    http://en.wikipedia.org/wiki/N-sphere#Generating_random_points
    """
    normal_deviates = norm(size=(N, dim))
    radius = np.sqrt((normal_deviates ** 2).sum(axis=0))
    points = normal_deviates / radius
    return points

# Number of hyperplanes
n = 10
# Dimension of space
d = 3

fig, ax = plt.subplots(subplot_kw=dict(projection='3d'))
points = points_on_sphere(n, d).T
uu, vv = np.meshgrid([-5, 5], [-5, 5], sparse=True)
colors = np.linspace(0, 1, len(points))
cmap = plt.get_cmap('jet')
for nhat, c in IT.izip(points, colors):
    u = (0, 1, 0) if np.allclose(nhat, (1, 0, 0)) else np.cross(nhat, (1, 0, 0))
    u /= math.sqrt((u ** 2).sum())
    v = np.cross(nhat, u)
    u = u[:, np.newaxis, np.newaxis]
    v = v[:, np.newaxis, np.newaxis]
    xx, yy, zz = u * uu + v * vv
    ax.plot_surface(xx, yy, zz, alpha=0.5, color=cmap(c))
ax.set_xlim3d([-5,5])
ax.set_ylim3d([-5,5])
ax.set_zlim3d([-5,5])        
plt.show()

或者,您可以使用 Till Hoffmann's pathpatch_2d_to_3d utility function 来避免繁琐的数学运算:

for nhat, c in IT.izip(points, colors):
    p = patches.Rectangle((-2.5, -2.5), 5, 5, color=cmap(c), alpha=0.5)
    ax.add_patch(p)
    pathpatch_2d_to_3d(p, z=0, normal=nhat)

ax.set_xlim3d([-5,5])
ax.set_ylim3d([-5,5])
ax.set_zlim3d([-5,5])        
plt.show()

【讨论】:

  • 我真的很喜欢你链接到 pathpatch_2d_to_3d ,谢谢!
【解决方案2】:

我建议你检查一下你的轴。您的计算使 Z 轴 方式 太大,这意味着您有一个荒谬的偏见观点。

首先检查您的法线是否均匀分布在圆上:

from __future__ import division
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

#Number of hyperplanes
n = 1000
#Dimension of space
d = 3

plt3d = plt.figure().gca(projection='3d')
for i in xrange(n):
    #Create random point on unit sphere
    v = np.random.normal(size = d)
    v = v/np.sqrt(np.sum(v**2))
    v *= 10

    plt3d.scatter(v[0], v[1], v[2])

plt3d.set_aspect(1)
plt3d.set_xlim(-10, 10)
plt3d.set_ylim(-10, 10)
plt3d.set_zlim(-10, 10)

plt.show()

然后检查您的飞机是否正确创建:

from __future__ import division
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

#Number of hyperplanes
n = 1
#Dimension of space
d = 3

plt3d = plt.figure().gca(projection='3d')
for i in xrange(n):
    #Create random point on unit sphere
    v = np.random.normal(size = d)
    v = v/np.sqrt(np.sum(v**2))
    v *= 10

    # create x,y
    xx, yy = np.meshgrid(np.arange(-5,5,0.3), np.arange(-5,5,0.3))
    xx = xx.flatten()
    yy = yy.flatten()
    z = (-v[0] * xx - v[1] * yy)/v[2]

    # Hack to keep the plane small
    filter = xx**2 + yy**2 + z**2 < 5**2
    xx = xx[filter]
    yy = yy[filter]
    z = z[filter]

    # plot the surface
    plt3d.scatter(xx, yy, z, alpha = 0.5)

    for i in np.arange(0.1, 1, 0.1):
        plt3d.scatter(i*v[0], i*v[1], i*v[2])

plt3d.set_aspect(1)
plt3d.set_xlim(-10, 10)
plt3d.set_ylim(-10, 10)
plt3d.set_zlim(-10, 10)

plt.show()

那么你就可以看到你其实已经得到了很好的结果!

from __future__ import division
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

#Number of hyperplanes
n = 100
#Dimension of space
d = 3

plt3d = plt.figure().gca(projection='3d')
for i in xrange(n):
    #Create random point on unit sphere
    v = np.random.normal(size = d)
    v = v/np.sqrt(np.sum(v**2))
    v *= 10

    # create x,y
    xx, yy = np.meshgrid(np.arange(-5,5,0.3), np.arange(-5,5,0.3))
    xx = xx.flatten()
    yy = yy.flatten()
    z = (-v[0] * xx - v[1] * yy)/v[2]

    # Hack to keep the plane small
    filter = xx**2 + yy**2 + z**2 < 5**2
    xx = xx[filter]
    yy = yy[filter]
    z = z[filter]

    # plot the surface
    plt3d.scatter(xx, yy, z, alpha = 0.5)

plt3d.set_aspect(1)
plt3d.set_xlim(-10, 10)
plt3d.set_ylim(-10, 10)
plt3d.set_zlim(-10, 10)

plt.show()

【讨论】:

    【解决方案3】:

    寻找不是一切。下次你最好测量一下:-]。它似乎是非随机分布的,因为您没有固定轴。因此,您会看到一架升到天空的主要飞机,而其余飞机由于规模,看起来非常相似,而不是随机分布。

    这段代码怎么样:

    from __future__ import division
    import numpy as np
    import matplotlib.pyplot as plt
    from mpl_toolkits.mplot3d import Axes3D
    
    #Number of hyperplanes
    n = 20
    #Dimension of space
    d = 3
    
    plt3d = plt.figure().gca(projection='3d')
    for i in xrange(n):
        #Create random point on unit sphere
        v = np.random.normal(size = d)
        v = v/np.sqrt(np.sum(v**2))
        # create x,y
        xx, yy = np.meshgrid(range(-1,1), range(-1,1))
        z = (-v[0] * xx - v[1] * yy)/v[2]
        # plot the surface
        plt3d.plot_surface(xx, yy, z, alpha = 0.5)
    
    plt3d.set_xlim3d([-1,1])
    plt3d.set_ylim3d([-1,1])
    plt3d.set_zlim3d([-1,1])
    plt.show()
    

    它并不完美,但现在看起来更加随机......

    【讨论】:

    • 我用来在球体上选择随机点的方法在 mathworld.wolfram.com/SpherePointPicking.html 的末尾。这然后将平面定义为math.stackexchange.com/questions/952525/… 显然我做错了什么。
    • 嗯,有道理,但说实话,我真的不明白它是如何成为均匀分布的......
    • @Jendas 哈哈!看起来你比我早一个小时就得到了答案,而不是 3 分钟!那么,由此产生的输出是a bit of a mess. 这也是你得到的吗?
    • 是的,我一直在编辑已删除的答案,终于到了这里。我也弄得一团糟,但我专注于尽可能保持原始代码不变,以便提问者更好地理解......另外,他质疑平面的均匀分布,这很重要。
    【解决方案4】:

    我试过这个,也许这是创建统一平面的更好方法。我对球坐标系随机取两个不同的角度,并将其转换为笛卡尔坐标,得到平面的法向量。此外,当您绘图时,您应该注意平面的中点不在原点上。

    import numpy as np
    import matplotlib.pyplot as plt
    from mpl_toolkits.mplot3d import Axes3D
    
    fig = plt.figure()
    ax = Axes3D(fig)
    
    for i in range(20):
        theta = 2*np.pi*np.random.uniform(-1,1)    
        psi = 2*np.pi*np.random.uniform(-1,1)
        normal = np.array([np.sin(theta)*np.cos(psi),np.sin(theta)*np.sin(psi),
                           np.cos(theta)])
        xx, yy = np.meshgrid(np.arange(-1,1), np.arange(-1,1))
        z = (-normal[0] * xx - normal[1] * yy)/normal[2]
        ax.plot_surface(xx, yy, z, alpha=0.5)    
    

    【讨论】:

      猜你喜欢
      • 2021-02-07
      • 2012-10-01
      • 1970-01-01
      • 2013-01-27
      • 1970-01-01
      • 1970-01-01
      • 2016-07-20
      • 2021-05-08
      • 1970-01-01
      相关资源
      最近更新 更多