python - 如何在给定纬度和经度列表的球体上绘制热图？

Question

I have a list of over 500 points, given in Latitude and Longitudes. These points represent craters, and I want to plot a heatmap of these craters. For example, I want an area with a lot of craters to be considered "hot" and fewer craters to be "cold". I have looked at KDE using SciPy, and also tried using ListSliceDensityPlot3D in Mathematica, but I have been unable to create a graph that is adequate.

I converted each point from latitude/longitude into Cartesian [x,y,z] coordinates, and plotted them on the surface of a sphere, but I don't know how I would take the list of points and calculate the density in a given area, and then plot it on a 3D surface.

The idea being that I obtain a plot something like this image of Ceres!

在此先感谢，如果需要，请提出问题，如果我最初没有发布足够的信息，对不起。

标签： python3dwolfram-mathematicaheatmap

Answer 1

这是一种蛮力方法，但它在一定程度上起作用。如果您使网格非常精细或有数千个陨石坑，这将是有问题的。如果 bin 尺寸足够小，则表面距离和 3D 距离之间没有太大差异，所以我采用后者，因为它更容易计算，但可能需要更改它。

代码如下：

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


def random_point( r=1 ):
    ct = 2*np.random.rand() - 1
    st = np.sqrt( 1 - ct**2 )
    phi = 2* np.pi *  np.random.rand()
    x = r * st * np.cos( phi)
    y = r * st * np.sin( phi)
    z = r * ct
    return np.array( [x, y, z ] )

def near( p, pntList, d0 ):
    cnt=0
    for pj in pntList:
        dist=np.linalg.norm( p - pj )
        if dist < d0:
            cnt += 1 - dist/d0
    return cnt


"""
https://stackoverflow.com/questions/22128909/plotting-the-temperature-distribution-on-a-sphere-with-python
"""

pointList = np.array([ random_point( 10.05 ) for i in range( 65 ) ] )

fig = plt.figure()
ax = fig.add_subplot( 1, 1, 1, projection='3d')

u = np.linspace( 0, 2 * np.pi, 120)
v = np.linspace( 0, np.pi, 60 )

# create the sphere surface
XX = 10 * np.outer( np.cos( u ), np.sin( v ) )
YY = 10 * np.outer( np.sin( u ), np.sin( v ) )
ZZ = 10 * np.outer( np.ones( np.size( u ) ), np.cos( v ) )

WW = XX.copy()
for i in range( len( XX ) ):
    for j in range( len( XX[0] ) ):
        x = XX[ i, j ]
        y = YY[ i, j ]
        z = ZZ[ i, j ]
        WW[ i, j ] = near(n p.array( [x, y, z ] ), pointList, 3)
WW = WW / np.amax( WW )
myheatmap = WW

# ~ ax.scatter( *zip( *pointList ), color='#dd00dd' )
ax.plot_surface( XX, YY,  ZZ, cstride=1, rstride=1, facecolors=cm.jet( myheatmap ) )
plt.show() 

结果是这样的：

你也可以修改距离函数来考虑陨石坑的大小，也许。