使用matplotlib创建全天投影时的奇怪（糟糕？）行为

Question

我正在尝试制作一个“全天空”密度图，这个图在RA（即0到360度）上是完整的，但在DEC上是不完整的（假设范围是从-45到90度）。如果我不使用任何投影直接绘制，这样是可以的，但当我尝试使用“莫尔维德投影”时，结果却没有恢复到我输入的数据。不过，如果我稍微修改一下代码，就能得到我预期的效果（不过，对于这个变化我没有一个合理的解释，正如你在例子中看到的）。

让我们来看一个完整的例子和它的输出，这样会更清楚：

import numpy as np
import matplotlib.pyplot as plt
import matplotlib.backends.backend_agg
from math import pi

#array between 0 and 360 deg
RA = np.random.random(10000)*360
#array between -45 and 90 degrees. By construction!
DEC= np.random.random(10000)*135-45

fig = plt.Figure((10, 4.5))
ax = fig.add_subplot(111,projection='mollweide')
ax.grid(True)
ax.set_xlabel('RA')
ax.set_ylabel('DEC')

ax.set_xticklabels(np.arange(30,331,30))
hist,xedges,yedges = np.histogram2d(DEC,RA,bins=[90,180],range=[[-90,90],[0,360]])
#TO RECOVER THE EXPECTED BEHAVIOUR, I HAVE TO CHANGE -90 FOR -80 IN THE PREVIOUS LINE:
#hist,xedges,yedges = np.histogram2d(DEC,RA,bins=[90,180],range=[[-80,90],[0,360]])
#I DO NOT WHY!

extent = (-pi,pi,-pi/2.,pi/2.)
image = ax.imshow(hist,extent=extent,clip_on=False,aspect=0.5,origin='lower')

cb = fig.colorbar(image, orientation='horizontal')
canvas = matplotlib.backends.backend_agg.FigureCanvasAgg(fig)

fig.canvas.print_figure("image1.png")

输出的图像是：

[由于我在这里是新手，所以不能发布图片，我会发一个链接，如果链接不行，请给我发邮件，我可以分享一个Dropbox文件夹，里面有图片；）]

我得到的输出图像

在这里你可以清楚地看到RA的范围是正确的，确实在0到360之间，但DEC的范围是-35到90，而不是-45到90。到目前为止，我不明白为什么少了10度。

不过，如果我稍微修改一下代码，把这一行：

hist,xedges,yedges = np.histogram2d(DEC,RA,bins=[90,180],range=[[-90,90],[0,360]]

改成：

hist,xedges,yedges = np.histogram2d(DEC,RA,bins=[90,180],range=[[-80,90],[0,360]]

我就得到了我认为应该得到的结果，这个图是：

输出图像2

[再次说明，如果链接不行，请告诉我，我可以和你分享一个Dropbox文件夹]

在这个图中，DEC的范围现在是从-45到90，正如我所创建的那样。

不过，把-90改成-80似乎没有道理（我觉得）。所以可能我在某些地方做错了，但我现在没注意到，或者我对代码有误解，或者matplotlib里有个奇怪的bug？？

任何帮助、提示或纠正都将非常感谢。

爱德华多

Answer 1

如果你不介意使用一个外部的工具包，你可以用healpy来实现，它提供了一种叫做Mollweide投影的方式，用于处理Healpix的天空像素化：

https://github.com/healpy/healpy

这里有一个和你的脚本类似的例子：

https://gist.github.com/1215159

关于healpix的更多信息：

http://healpix.jpl.nasa.gov/html/intro.htm

输出的图片：

Answer 2

如果这对其他人有帮助的话，这是我代码的“修正版本”，它的输出结果是这个 图片。主要的改动是使用了 pcolormesh 而不是 imshow（正如 @Joe 建议的那样）：

import numpy as np
import matplotlib.pyplot as plt
import matplotlib.backends.backend_agg

#array between 0 and 360 deg
#CAVEAT: it seems that is needed an array from -180 to 180, so is just a
#shift in the coordinates
RA = np.random.random(10000)*360-180
#array between -45 and 90 degrees
DEC= np.random.random(10000)*135-45

fig = plt.Figure((10, 5))
ax = fig.add_subplot(111,projection='mollweide')

ax.set_xlabel('RA')
ax.set_ylabel('DEC')
ax.set_xticklabels(np.arange(30,331,30))

hist,xedges,yedges = np.histogram2d(DEC,RA,bins=[60,40],range=[[-90,90],[-180,180]])

X,Y = np.meshgrid(np.radians(yedges),np.radians(xedges))

image = ax.pcolormesh(X,Y,hist)
ax.grid(True)
cb = fig.colorbar(image, orientation='horizontal')
canvas = matplotlib.backends.backend_agg.FigureCanvasAgg(fig)
fig.canvas.print_figure("image4.png")