천리안 위성 2A호(GK2A)의 데이터 표출 [3]

천리안 2A 위성의 영상 데이터를 표출하는 방법과 예제를 소개하는 게시물들을 제가 지난 해(2020)에 2회에 걸쳐서 게재했던 바 있습니다.

 

천리안 위성 2A호(GK2A)의 데이터 표출 [1]

천리안 위성 2A호(GK2A)의 데이터 표출 [2]

 

여기서는 적외선 영상인 IR 11.2 um 데이터를 담고 있는 netCDF 포맷의 파일을 사용하였습니다. 그리고 영상 데이터의 영역 커버리지가 전구(Global), 동아시아, 한반도 등 세가지 종류로 제공되는데, 동아시아 및 한반도 영역 데이터를 대상으로 하였고 전구 맵(Global Map) 데이터의 경우는 보류했었는데요. 오늘은 이 전구 맵 데이터를 대상으로 표출 방법 및 예제를 소개해보고자 합니다. 이번에 예제로 사용할 데이터 파일의 이름은 다음과 같습니다.

 

gk2a_ami_le1b_ir112_fd020ge_202002280300.nc

 

파일 이름에서 알 수 있듯이 2020년 2월 28일 03시 00분 UT에 IR 11.2 um 파장으로 얻어진 전구(Global) 영역에 대한 데이터를 담은 netCDF 형식의 파일입니다. 다만 이 파일은 천리안 2A 자료검색 시스템 웹페이지에 가입을 하여 요청을 해서 받은 것이기 때문에, 제가 여기서 이 파일을 배포하는 것은 적절하지 않을 것 같습니다. 따라서 데이터 파일이 필요하신 분들은 웹페이지에 직접 가입을 해서 받으시길 권장드립니다. 그리고 제가 편의상 여기서는 IR 11.2 um 데이터로 한정하였지만, 이번에 소개될 내용은 다른 종류의 데이터에 대해서도 거의 비슷하게 적용될 수 있음을 미리 언급해둡니다.

 

그리고 한가지 더 유념해야 할 것은, 나중에 언급하겠지만 천리안 2A 전구 영역 데이터의 경우 GOES-R 투영법을 기반으로 한다는 것입니다. 정지궤도 위성들의 일반적인 특성이죠. 따라서 제가 지금부터 소개할 표출 과정에서는 바탕 지도를 당연히 GOES-R 투영법으로 설정해야 하는데요. IDL에서 GOES-R 투영법이 지원되기 시작한 것이 IDL 8.6.1 버전부터입니다. 따라서 오늘 소개할 내용은 IDL 8.6.1 이전의 버전에서는 적용되지 않는다는 점을 미리 언급해둡니다. 그러면 본격적으로 시작해보겠습니다. 먼저 다음과 같이 파일명을 문자값으로서 변수에 담아둡시다.

file = 'gk2a_ami_le1b_ir112_fd020ge_202002280300.nc'

 

그러면 우선 이 파일 내에 포함되어 있는 변수(Variable)들의 목록 정보부터 먼저 살펴보겠습니다. 이를 위하여 다음과 같이 NCDF_LIST 명령을 사용합니다.

NCDF_LIST, file, /VARIABLES

이 명령에 의하여 출력된 내용은 원래 좀 길지만 앞부분만 보면 다음과 같습니다.

 

 gk2a_ami_le1b_ir112_fd020ge_202002280300.nc
     # dimensions: 8
     # Variables: 65
     # Global attributes: 53
     The unlimited dimension is 2
 
     Variables
           0  image_pixel_values:                        UINT(5500,5500) = UINT(dim_image_x,dim_image_y)

 

이 내용에 의하면 변수는 65개이고 그 중 첫번째 변수가 바로 5500x5500의 구조를 갖는 2차원 배열 데이터임을 알 수 있습니다. 이것이 결국 이미지 데이터에 해당되며 그 명칭은 'image_pixel_values'입니다. 따라서 나중에 이 배열 데이터를 추출할 때 이 명칭을 이용할 것입니다. 그리고 전역속성(Global Attribute)이 53개의 항목들로 구성되어 있음을 확인할 수 있습니다. 그래서 이 전역속성들의 세부 내용을 살펴보기 위하여 다음과 같이 NCDF_LIST 명령을 한번 더 사용합니다.

 

NCDF_LIST, file, /GATT

 

이 명령에 의하여 출력된 내용은 다음과 같습니다.

 

     Global Attributes
           0  satellite_name: GK-2A
           1  instrument_name: AMI
           2  data_processing_center: NMSC
           3  data_processing_mode: operation
           4  channel_spatial_resolution: 2.0
           5  channel_center_wavelength: 11.2
           6  scene_acquisition_time: 20200228_030032
           7  mission_reference_time: 20200228_030000
           8  file_generation_time: 20200228_030956
           9  file_name: gk2a_ami_le1b_ir112_fd020ge_202002280300.nc
          10  file_format_version: 1.0.0_20181120
          11  geometric_correction_sw_version: GK2_INRSM_V1.3
          12  star_catalog_version: 1.0.0_20170303
          13  landmark_catalog_version: 1.0.0_20190502
          14  calibration_table_version: v.3.0_20190415
          15  DN_to_Radiance_Gain: -0.021674486
          16  DN_to_Radiance_Offset: 176.71344
          17  Teff_to_Tbb_c0: -0.24911172
          18  Teff_to_Tbb_c1: 1.0012117
          19  Teff_to_Tbb_c2: -1.1316796e-06
          20  light_speed: 2.9979246e+08
          21  Plank_constant_h: 6.6260696e-34
          22  Boltzmann_constant_k: 1.3806488e-23
          23  total_pixel_data_size: 60500000
          24  number_of_total_swaths: 23
          25  number_of_columns: 5500
          26  number_of_lines: 5500
          27  observation_mode: FD
          28  observation_start_time: 6.3613083e+08
          29  observation_end_time: 6.3613138e+08
          30  time_synchro_obt: 6.3613083e+08
          31  time_synchro_utc: 6.3613083e+08
          32  projection_type: GEOS
          33  sub_longitude: 2.2375121
          34  cfac: 20425339.
          35  lfac: -20425339.
          36  coff: 2750.5000
          37  loff: 2750.5000
          38  nominal_satellite_height: 42164000.
          39  earth_equatorial_radius: 6378137.0
          40  earth_polar_radius: 6356752.3
          41  earth_ellipsoid_center: 0.0000000
          41  earth_ellipsoid_center: 0.0000000
          41  earth_ellipsoid_center: 0.0000000
          42  image_upperleft_latitude: -1.0000000e+14
          43  image_upperleft_longitude: -1.0000000e+14
          44  image_lowerright_latitude: -1.0000000e+14
          45  image_lowerright_longitude: -1.0000000e+14
          46  image_center_latitude: 0.0000000
          47  image_center_longitude: 2.2375121
          48  image_upperleft_x: -0.15397200
          49  image_upperleft_y: 0.15397200
          50  image_lowerright_x: 0.15397200
          51  image_lowerright_y: -0.15397200
          52  resampling_kernel_type: SINC

 

이와 같이 총 53종의 정보들이 출력되는데, 이들 중에는 표출 작업에서 필요한 것들이 여럿 있습니다. 가장 눈에 띄는 것은 projection type 정보인데 값을 보면 'GEOS'라고 되어 있습니다. 이것은 Geostationary 투영법이라는 의미인데 사실 이 GEOS 투영법은 GOES-R 투영법과는 엄연히 다르다는 것을 유의해야 합니다. 일단 Geostationary 투영법에 관한 자세한 내용은 아래 링크를 참조하시면 됩니다.

 

https://proj.org/operations/projections/geos.html

 

GEOS(Geostationary) 투영법과 GOES-R 투영법 사이의 주요 차이점은 두가지입니다. 첫번째 차이점은 지상고도의 값인데 GEOS의 경우는 35785831m이고 GOES-R의 경우는 35786023m입니다. 수치상으로 다른 것은 맞지만 그 차이는 사실 미미한 수준이기도 합니다. 두번째 차이점은 스캔 방향의 순서인데 세부 내용은 약간 복잡하지만 결과적으로는 이미지 데이터의 X축 및 Y축 방향이 서로 반대가 됩니다. 그래서 나중에 나오겠지만 이미지 배열을 추출한 후에 REVERSE 함수를 사용하여 X 및 Y 방향을 반대로 뒤집는 과정이 추가될 것입니다. 따라서 GEOS 투영법과 GOES-R 투영법은 원론적으로는 당연히 서로 차이가 있지만, IDL에서는 GOES-R 투영법을 사용하면서 이미지 배열을 뒤집어주면 GEOS 투영법과 거의 비슷한 효과를 얻을 수 있다고 봐도 될 것 같습니다. 따라서 여기서는 바탕 지도를 표출하는데 있어서 GOES-R 투영법을 사용하기로 하겠습니다.

 

그 외에도 바탕 지도의 투영법과 관련된 여러가지 중요한 세부 정보들이 있습니다. 그래서 이후의 작업에서 중요하게 사용될 변수 및 속성 항목들을 골라서 다음과 같이 NCDF 관련 명령들을 사용하여 배열 또는 변수의 형태로 가져옵시다. 각 항목에 대해서는 이후의 표출과정에서 자세히 살펴보겠습니다.

 

id = NCDF_OPEN(file)
NCDF_VARGET, id, 'image_pixel_values', img
HELP, img
NCDF_ATTGET, id, 'image_center_longitude', value, /GLOBAL
clon = value*!radeg
NCDF_ATTGET, id, 'image_center_latitude', value, /GLOBAL
clat = value*!radeg
NCDF_ATTGET, id, 'nominal_satellite_height', hgt, /GLOBAL
PRINT, clon, clat, hgt
NCDF_ATTGET, id, 'earth_equatorial_radius', radius_eq, /GLOBAL
NCDF_ATTGET, id, 'earth_polar_radius', radius_pol, /GLOBAL
PRINT, radius_eq, radius_pol
NCDF_ATTGET, id, 'image_lowerright_x', x_lr, /GLOBAL
PRINT, x_lr
NCDF_ATTGET, id, 'image_upperleft_y', y_lr, /GLOBAL
PRINT, y_lr
NCDF_CLOSE, id

 

이와 같이 꽤 많은 항목의 정보들을 추출하면서 각 항목의 내용도 출력하도록 하였습니다. 그런데 이렇게 추출한 정보들 중에서 실제 표출에서 필요한 항목들은 img, clon, x_lr, y_lr 정도 뿐입니다. 그 외의 것들은 정보를 확인하는 차원에서 추출해본 것이고 실제 표출 과정에서는 사용되지 않습니다. 이 부분은 표출 과정을 진행하면서 살펴보기로 하겠습니다. 먼저 바탕이 될 GOES-R 투영법의 지도를 다음과 같이 표출해봅시다.

 

win = WINDOW(DIMENSIONS=[800, 800], /NO_TOOLBAR)
m = MAP('GOES-R', CENTER_LONGITUDE=clon, $
  HORIZON_LINESTYLE=6, MARGIN=0.05, /CURRENT)

 

이와 같이 MAP 함수를 사용하면서 투영법을 'GOES-R'로 명시하였습니다. 이 과정에서 따로 설정해야 할 중요한 속성은 바로 중심 경도입니다. IDL에서 GOES-R 투영법을 사용하는 방법에 관해서는 제가 얼마전에 관련 게시물을 올린 바 있고 거기서도 언급을 했듯이, 이 투영법은 정지궤도 위성이 바라보는 지구 타원체의 모습을 표출하는 목적이고 위성의 고도 및 지구 타원체에 대한 기하학적 수치들에 대한 디폴트 값들이 정의되어 있습니다. 여기서는 이 디폴트 값들을 그대로 사용하면 되고, 다만 중심 경도만 따로 설정해주면 됩니다. 이를 위하여 앞서 clon이라는 변수에 중심 경도 값(128.2도)을 담아둔 것입니다.

 

바탕 지도를 표출한 다음 과정은 이미지 데이터를 중첩하는 것입니다. 이미지 데이터는 앞서 img라는 배열로 추출했는데 일단 다음과 같이 배열 img에 대한 기본적인 정보를 출력해봅시다.

 

HELP, img

PRINT, MIN(img), MAX(img)

 

출력된 내용은 다음과 같습니다.

 

IMG             UINT      = Array[5500, 5500]
     110   32768

 

이 내용을 보면 img는 5500x5500의 구조인 2차원 배열이고 최소 및 최대 값은 110과 32768로 나옵니다. 그런데 제가 보기에는 아마도 32768은 의미가 없는 값들인 것으로 보여서 다음과 같이 이러한 값들에 대해서는 NaN으로 처리를 해보았습니다(물론 제 짐작이 틀릴 수도 있습니다). 그리고 이와 더불어 REVERSE 함수를 사용하여 배열을 뒤집는 처리도 함께 하였는데, 원래 이미지 데이터를 그대로 추출한 상태는 우리가 보기에 Y축에 대하여 뒤집어진 상태이기 때문입니다.

 

img = REVERSE(FLOAT(img), 2)
ww = WHERE(img GT 30000, count)
IF count NE 0 THEN img[ww] = !values.f_nan
PRINT, MIN(img), MAX(img)

 

이러한 처리를 거친 후 다시 img의 최소 및 최대 값을 확인해보면 다음과 같습니다.

 

      110.000      7972.00

 

이제 IMAGE 함수를 사용하여 배열 img를 이미지의 형태로 중첩해야 합니다. 그 과정은 다음과 같이 처리하였습니다.

 

im = IMAGE(img, IMAGE_LOCATION=[-x_lr, -y_lr], $
  IMAGE_DIMENSION=[x_lr*2, y_lr*2], GRID_UNITS=1, $

  RGB_TABLE=0, /OVERPLOT)

 

여기서는 IMAGE_LOCATION 및 IMAGE_DIMENSION 속성을 사용하여 이미지의 위치 및 크기를 지정하였는데, 이 과정에서는 x_lr 및 y_lr이 사용되었습니다. 이 변수들은 앞서 NCDF 파일에서 추출한 image_lowerright_x 및 image_upperleft_y라는 전역 속성 항목에 해당됩니다. 쉽게 생각하면 이미지의 X 방향 및 Y 방향 끝 지점의 좌표값이라고 보면 됩니다. 그리고 이러한 값들은 경위도가 아닌 거리에 해당되는 값들이기 때문에 GRID_UNITS 속성의 값을 1로 설정한 것입니다. 다만 여기서 주목해볼 것은 바로 x_lr, y_lr의 값들인데 앞서 이미 출력을 해봤듯이 그 값은 둘 다 0.153972라고 나옵니다. 이 값은 실제로 scan-angle radian에 해당됩니다. 그리고 IDL의 MAP 함수로 표출된 바탕 지도의 투영법이 GOES-R인 상태에서는 이미지 데이터를 중첩할 때 그 위치 및 크기에 해당되는 좌표값을 이러한 단위로 받도록 되어 있습니다. 이 내용은 IDL 도움말의 MAP 함수에 관한 섹션에서 찾아볼 수 있는데 해당 부분은 다음과 같습니다.

 

 

이러한 이유로 인하여 위와 같이 이미지의 위치 및 크기를 설정하였다고 보면 됩니다. 이제 이미지 중첩까지 되었으므로 남은 과정은 대륙경계선 및 경위도 격자선을 표시하는 것입니다. 그 과정은 다음과 같이 처리해보았습니다.

 

mc = MAPCONTINENTS(/HIRES, COLOR='orange red')
m.MapGrid.Linestyle = 0
m.MapGrid.COLOR = 'deep sky blue'
m.MapGrid.GRID_LONGITUDE = 15
m.MapGrid.GRID_LATITUDE = 15

 

이제 다 되었습니다. 표출된 그림을 보면 다음과 같습니다.

 

 

그림을 보면 바탕 지도와 이미지가 잘 중첩이 된 것 같습니다. 따라서 다른 종류의 전구 맵 데이터들의 경우도 유사한 방식으로 처리하여 표출하면 될 것으로 생각됩니다. 그리고 아마도 천리안 외 다른 정지궤도 위성들의 전구 맵 데이터의 경우도 세부적으로 약간 차이가 날 수도 있겠지만 전반적인 처리 및 표출의 기본적인 과정은 이와 같다고 보면 되지 않을까 싶습니다.