將旋轉應用到healpy中的HEALPix地圖

我有一個HEALPix地圖，我已經閱讀過使用healpy，但它在銀河座標系中，我需要它在天體/赤道座標系中。有誰知道轉換地圖的簡單方法嗎？將旋轉應用到healpy中的HEALPix地圖

我已經嘗試使用healpy.Rotator從（l，b）轉換爲（phi，theta），然後使用healpy.ang2pix重新排列像素，但該圖仍然看起來很奇怪。

如果有類似於Rotator的功能，您可以這樣稱呼它：map = AnotherRotator(map,coord=['G','C'])。任何人都知道任何這樣的功能？

感謝，

亞歷

2014-07-08 aim

的可能的複製[如何轉換和保存healpy映射到不同的座標系？（https://stackoverflow.com/questions/44443498/how-to-convert-and- save-healpy-map-to-different-coordinate-system） –

我意識到，這是問很久以前的事，但我有同樣的問題我自己這一週，發現您的文章。我發現了一些潛在的解決方案，所以我會分享一些其他人遇到此問題並發現它有用。

解決方案1 ：這種取決於您的數據進入的格式。礦進入（theta，phi）網格。

import numpy as np 
import healpy as H 
map = <your original map> 
nside = <your map resolution, mine=256> 
npix = H.nside2npix(nside) 
pix = N.arange(npix) 
t,p = H.pix2ang(nside,pix) #theta, phi 

r = H.Rotator(deg=True, rot=[<THETA ROTATION>, <PHI ROTATION>]) 

map_rot = np.zeros(npix) 

for i in pix: 
    trot, prot = r(t[i],p[i]) 
    tpix = int(trot*180./np.pi) #my data came in a theta, phi grid -- this finds its location there 
    ppix = int(prot*180./np.pi) 
    map_rot[i] = map[ppix,tpix] #this being the rright way round may need double-checking

解決方案2：還沒有完全完成測試這一點，但只是通過它來了以後做上述煩人的工作...

map_rot = H.mollview(map,deg=True,rot=[<THETA>,<PHI>], return_projected_map=True)

這給2D numpy的陣列。我很想知道如何將其轉換回healpix映射...

來源

2014-10-31 18:37:04

我發現了另一個潛在的解決方案，經過幾個月搜索和搜索。我還沒有測試過，所以請小心！

掃羅的解決方案2，上面是關鍵（偉大的建議！）

基本上，你是在reproject_to_healpix功能相結合的healpy.mollview功能（gnomview，cartview和orthview工作，以及） reproject包（http://reproject.readthedocs.org/en/stable/）。

由此產生的地圖適用於我的角度尺度，但我不能說與其他方法相比轉換有多精確。

-----基本輪廓----------

步驟1：讀入的地圖，並經由cartview使矩形陣列。正如掃羅上面所指出的，這也是進行輪換的一種方式。如果您只是進行標準的旋轉/座標轉換，那麼您只需要輸入coord關鍵字。從富到銀河座標，設定coord = ['C','G']

map_Gal = hp.cartview(map_Cel, coord=['C','G'], return_projected_map=True, xsize=desired_xsize, norm='hist',nest=False)

步驟2：寫模板全天空FITS報頭（如在下面的例子）。我寫我的像我所期望的HEALPix地圖一樣具有平均像素級。

步驟3：使用reproject.transform_to_healpix

reproject包括用於映射「正常」數組的函數（或FITS文件）插入HEALPix投影。將它與返回由healpy.mollview/cartview/orthview/gnomview創建的數組的能力結合起來，並且可以將一個座標系（天體）的HEALPix映射旋轉到另一個座標系（星系）。

map_Gal_HP, footprint_Gal_HP = rp.reproject_to_healpix((map_Gal, target_header), coord_system_out= 'GALACTIC', nside=nside, nested=False)

它實際上歸結爲這兩個命令。但是你必須製作一個模板標題，給出與你想要製作的中間全天空地圖相對應的像素比例和大小。

-----全部工作例子（IPython的筆記本格式+ FITS樣本數據）------

https://github.com/aaroncnb/healpix_coordtrans_example/tree/master

的代碼應該運行得非常快，但是那是因爲地圖嚴重退化。我爲我的NSIDE 1024和2048地圖做了同樣的工作，花了大約一個小時就花費了。

------之前和之後的圖像------

來源

2016-04-25 05:26:30

對不起，我盲目地做了一個與問題相反的例子：我的例子將天體轉換爲銀河座標。請求了天體銀河。但是，您可以輕鬆修改命令以將星系轉換爲天體。 cartview'中的coord = ['C'，'G']''關鍵字變爲''coord = ['G'，'C']''。示例模板標題的CTYPE關鍵字將被修改爲「RA --- CAR」和「DEC-CAR」，而「COORDSYS」關鍵字將變爲「'icrs' '' –

這個功能似乎這樣的伎倆（相當緩慢的，但應該是比for循環）好：

def rotate_map(hmap, rot_theta, rot_phi): 
    """ 
    Take hmap (a healpix map array) and return another healpix map array 
    which is ordered such that it has been rotated in (theta, phi) by the 
    amounts given. 
    """ 
    nside = hp.npix2nside(len(hmap)) 

    # Get theta, phi for non-rotated map 
    t,p = hp.pix2ang(nside, np.arange(hp.nside2npix(nside))) #theta, phi 

    # Define a rotator 
    r = hp.Rotator(deg=False, rot=[rot_phi,rot_theta]) 

    # Get theta, phi under rotated co-ordinates 
    trot, prot = r(t,p) 

    # Interpolate map onto these co-ordinates 
    rot_map = hp.get_interp_val(hmap, trot, prot) 

    return rot_map

數據從PyGSM使用此給出以下：

hp.mollview(np.log(rotate_map(gsm.generated_map_data, 0,0)))