Переключение между координатами пикселя и физическими координатами в n-мерном изображении FITS с помощью Astropy - PullRequest
0 голосов
/ 28 марта 2019

У меня есть изображение FITS с заголовком, содержащим информацию о координатах вдоль каждой оси

NAXIS   =                    3
NAXIS1  =                  259
NAXIS2  =                  272
NAXIS3  =                  100
CDELT1  = -0.08333333330000001
CDELT2  =  0.08333333330000001
CDELT3  =                  0.1
CRPIX1  =                130.5
CRPIX2  =    132.1906015909634
CRPIX3  =                    0
CRVAL1  =                255.0
CRVAL2  =                 60.0
CRVAL3  =                15.45

Есть ли простой способ (например, функция Astropy), чтобы получить фактические (физические) координаты из пикселякоординаты?

И наоборот, есть ли в Astropy встроенная функция для получения координат ближайшего пикселя из фактических координат?

Редактировать: Я нашел функцию WCS.all_pix2world, но я не понимаю, как ее использовать, и я не нашел никаких примеров использования.

Из документации:

There are two accepted forms for the positional arguments:

        - 2 arguments: An N x naxis array of coordinates, and an origin.
        - more than 2 arguments: An array for each axis, followed by an origin. These arrays must be broadcastable to one another.

В первом случае, что такое N?Не будет ли origin всегда [0,0,0]?

1 Ответ

2 голосов
/ 29 марта 2019
  • N - это количество точек, для которых вы хотите выполнить преобразование, т.е. вы можете вычислить пиксельные или мировые координаты для массива из N точек.

  • origin просто вопрос соглашения о координате верхнего левого угла.В массивах Python используется индексирование на основе 0, поэтому оно равно 0. См. Пример ниже:

In [1]: from astropy.wcs import WCS                                                            

In [2]: from astropy.io import fits                                                            

In [3]: hdr = fits.Header.fromstring("""\ 
   ...: NAXIS   =                    3 
   ...: NAXIS1  =                  259 
   ...: NAXIS2  =                  272 
   ...: NAXIS3  =                  100 
   ...: CDELT1  = -0.08333333330000001 
   ...: CDELT2  =  0.08333333330000001 
   ...: CDELT3  =                  0.1 
   ...: CRPIX1  =                130.5 
   ...: CRPIX2  =    132.1906015909634 
   ...: CRPIX3  =                    0 
   ...: CRVAL1  =                255.0 
   ...: CRVAL2  =                 60.0 
   ...: CRVAL3  =                15.45 
   ...: """, sep='\n')                                                                         

In [4]: wcs = WCS(hdr)                                                                         

In [5]: wcs                                                                                    
Out[5]: 
WCS Keywords

Number of WCS axes: 3
CTYPE : ''  ''  ''  
CRVAL : 255.0  60.0  15.45  
CRPIX : 130.5  132.1906015909634  0.0  
PC1_1 PC1_2 PC1_3  : 1.0  0.0  0.0  
PC2_1 PC2_2 PC2_3  : 0.0  1.0  0.0  
PC3_1 PC3_2 PC3_3  : 0.0  0.0  1.0  
CDELT : -0.0833333333  0.0833333333  0.1  
NAXIS : 259  272  100

Теперь вы можете вычислить мировую координату верхнего левого угла:

In [6]: wcs.all_pix2world([[0, 0, 0]], 0)                                                      
Out[6]: array([[265.79166666,  49.06744987,  15.55      ]])

Здесь вы можете использовать массив индексов Nx3.

Так что, если у вас есть индексы пикселей и вы хотите преобразовать их в координаты неба, вам нужно использовать origin = 0, и то же самое для противоположного,преобразовать координаты неба в индексы пикселей с помощью wcs.all_world2pix.

Использование origin = 1 иногда может быть полезно, если у вас есть индексы пикселей, сохраненные в каталоге с использованием соглашения FITS / Fortran.

...