Любопытное (плохое?) Поведение, создающее проекции всего неба с помощью matplotlib - PullRequest
1 голос
/ 09 сентября 2011

Я пытаюсь создать график плотности "все небо", который является полным в RA (то есть от 0 до 360 градусов), но неполным в DEC (скажем, от -45 до 90 градусов). Если я построю это без какой-либо проекции, это нормально, но когда я пытаюсь построить с использованием проекции 'mollweide', я не восстанавливаю ввод, но если я делаю небольшое изменение в коде, я восстанавливаю ожидаемое поведение (однако я нет четкого объяснения этому изменению, как вы увидите в примере).

Давайте посмотрим на отдельный пример с более ясными выводами:

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], диапазон = [[ -90 , 90], [0,360]]

для

hist, xedges, yedges = np.histogram2d (DEC, RA, bins = [90,180], диапазон = [[ -80 , 90], [0,360]]

Я получаю то, что, как мне кажется, я должен получить, вот этот график:

Выходное изображение 2

[Опять же, если ссылка не работает, дайте мне знать, и я могу поделиться с вами папкой Dropbox]

где DEC теперь варьируется от -45 до 90, как и ожидалось, потому что я создал DEC таким образом.

Однако изменение -90 на -80 не имеет смысла (я думаю). Поэтому, вероятно, я делаю что-то не так, что сейчас не могу заметить, или я что-то неправильно понимаю в коде, или есть странная ошибка в matplotlib ??

Пожалуйста, любая помощь / подсказка / исправление будет очень признателен

Эдуардо

Ответы [ 2 ]

3 голосов
/ 14 сентября 2011

если вы не возражаете против зависимости от внешнего пакета, вы можете сделать это с помощью healpy, который обеспечивает проекцию Mollweide для пикселизации неба Healpix:

https://github.com/healpy/healpy

См.Пример, подобный вашему сценарию:

https://gist.github.com/1215159

Подробнее о healpix:

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

Выходное изображение: Example hitmap in Mollview

2 голосов
/ 14 сентября 2011

Если это полезно для кого-то другого, это «исправленная версия» моего кода, которая выдает на выходе image .Основным изменением является использование 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")
Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...