GOES 16/17 데이터의 표출

GOES 16 및 17은 미국 NASA/NOAA에서 운용중인 기상관측용 정지궤도 인공위성들입니다. GOES 16은 GOES-R 또는 GOES East라고도 부르며 경도 위치가 75.2W이고, GOES-17은 GOES-S 또는 GOES West라고도 부르며 경도 위치는 137.2W입니다. 따라서 아메리카 대륙 및 대평양 지역에 대한 기상 연구에 있어서 GOES 16 및 17의 데이터가 자주 사용됩니다. 오늘은 이 데이터 파일을 읽고 표출하는 예제를 살펴보고자 합니다. 어차피 데이터의 종류가 매우 많지만 그 중에서 GOES 17의 ABI(Advanced Baseline Image)에서 관측되는 복사량(Radiance) 데이터를 예제 데이터로 사용해보겠습니다. GOES 16/17의 데이터를 받을 수 있는 경로들은 몇가지가 있는데, 그 중에서 아마존(Amazon) 데이터 서비스를 이용하는 방법이 있습니다. 링크는 아래와 같습니다.

 

GOES-16/17 on Amazon Download Page (utah.edu)

 

이 웹페이지에서 각종 데이터 파일들을 비교적 편리하게 받을 수가 있습니다. 여기서 제가 받아본 파일의 이름은 다음과 같습니다.

 

OR_ABI-L1b-RadF-M6C13_G17_s20220150350320_e20220150359397_c20220150359437.nc

 

GOES 17의 ABI L1B Radiance 데이터는 총 16개의 밴드들로 구성되어 있는데 그 중 13번 밴드(10.35 um)에 해당되는 전지구(Full Disk) 데이터를 netCDF 형식의 파일로 받은 것입니다. 실제로 제가 이 파일을 받기 위하여 웹페이지에서 설정한 사항들은 다음 스크린샷에서 보실 수 있습니다. 날짜 및 시각은 UTC로 2022년 1월 15일 03시 50분입니다.

 

 

그러면 이 파일을 입수한 상태에서 이제 IDL의 netCDF 파일 읽기 기능을 사용하여 파일을 읽어보겠습니다. 일단 이 파일에 수록된 변수들(Variables)의 개요를 살펴보기 위하여 다음과 같이 NCDF_LIST 명령을 사용해봅시다.

 

file = 'OR_ABI-L1b-RadF-M6C13_G17_s20220150350320_e20220150359397_c20220150359437.nc'

NCDF_LIST, file, /VARIABLES, /VATT

 

이 내용을 실행하면 파일 내에 수록된 모든 변수들의 개요 및 세부 정보들이 대량으로 출력되는데, 이 내용을 통하여 우리가 데이터의 표출을 위하여 어떤 변수들을 읽어야 할 것인지를 대략 가늠할 수 있습니다. 가장 먼저 이미지 데이터에 해당되는 'Rad'라는 변수를 읽어야 합니다. 또한 이 변수에 종속된 하위속성(Attribute)에 해당되는 항목들(scale_factor 및 add_offset)도 함께 읽어야 합니다. 그 과정은 다음과 같습니다.

 

id = NCDF_OPEN(file)
NCDF_VARGET, id, 'Rad', img
HELP, img
PRINT, MIN(img), MAX(img)
NCDF_ATTGET, id, 'Rad', 'scale_factor', d_scl
NCDF_ATTGET, id, 'Rad', 'add_offset', d_off

 

이와 같이 변수 'Rad'에 해당되는 항목을 img라는 배열로 읽어서 HELP 및 PRINT로 관련 정보를 출력해보면 다음과 같습니다.

 

IMG             INT       = Array[5424, 5424]
     175    4095

 

이와 같이 5424x5424의 구조를 갖는 2차원 정수 배열이고, 최소 및 최대 값은 각각 175 및 4095라고 나옵니다. 사실 이렇게 읽어들인 img에 대해서는 후속 처리가 좀 더 필요합니다. 왜냐하면 앞서 NCDF_LIST 명령으로 출력된 정보들 중에서 변수 'Rad'에 관한 부분에 매우 중요한 정보가 있기 때문입니다.

 

           0  Rad:                                              INT(5424,5424) = INT(x,y)
                 0  _FillValue:     4095
                 1  long_name: ABI L1b Radiances
                 2  standard_name: toa_outgoing_radiance_per_unit_wavenumber
                 3  _Unsigned: true
                 4  sensor_band_bit_depth: 
                 5  valid_range:        0
                 5  valid_range:     4094
                 6  scale_factor:     0.0457289
                 7  add_offset:      -1.64430
                 8  units: mW m-2 sr-1 (cm-1)-1
                 9  resolution: y: 0.000056 rad x: 0.000056 rad

 

즉 이 내용에 의하면 유효한 값의 범위는 0~4094이며 4095는 의미없는 값이 됩니다. 그리고 유효한 정수형 값들에 대해서는 scale_factor에 해당되는 값을 곱하고 add_offset에 해당되는 값을 더하는 변환까지 적용해줘야 실수형의 실제 데이터 값이 됩니다. 따라서 의미없는 값(4095)을 필터링하고 유의미한 값들에 대한 변환식 적용이 필요한데 그 과정은 다음과 같습니다.

 

ww = WHERE(img LT 0 OR img GT 4094, cbad)
PRINT, 'number of bad values :', cbad
img = FLOAT(img)
IF cbad NE 0 THEN img[ww] = !values.f_nan
img = img*d_scl+d_off
HELP, img
PRINT, MIN(img), MAX(img)

 

여기서는 먼저 유효 범위를 벗어나는 값들을 WHERE 함수로 탐색하여 그러한 의미없는 값들을 NaN으로 대체합니다. 그 다음에는 scale 및 offset 값을 적용하여 실수형 데이터 값으로 변환해줍니다. 이러한 처리를 거친 img에 대한 정보를 HELP 및 PRINT로 출력해본 결과는 다음과 같습니다.

 

IMG             FLOAT     = Array[5424, 5424]
      6.35826      113.638

 

이제 img는 실수형의 실제 데이터 값들로 구성된 배열이 되었으며, 이 배열 내에서는 의미없는 값들이 NaN으로 존재하게 됩니다. 좀 번거롭게 보일 수도 있겠지만 이러한 처리 과정은 반드시 필요합니다. 그리고 netCDF 파일로부터 읽어야 할 정보들은 아직 남아있습니다. 바로 변수 'x' 및 'y'인데요. 읽어들이는 과정은 다음과 같습니다.

 

NCDF_VARGET, id, 'x', x
NCDF_ATTGET, id, 'x', 'scale_factor', x_scl
NCDF_ATTGET, id, 'x', 'add_offset', x_off
x = x*x_scl+x_off
PRINT, MIN(x), MAX(x)
NCDF_VARGET, id, 'y', y
NCDF_ATTGET, id, 'y', 'scale_factor', y_scl
NCDF_ATTGET, id, 'y', 'add_offset', y_off
y = y*y_scl+y_off
PRINT, MIN(y), MAX(y)

 

여기서는 이미지를 구성하는 격자들의 X 및 Y 방향 위치값들(각각 5424개)을 읽고, 각각에 대하여 scale 및 offset을 적용하여 변환해주는 과정까지 처리하였습니다. X 및 Y 방향 모두 격자값들의 최소 및 최대 값은 다음과 같습니다.

 

    -0.151844     0.151844

 

이 값들의 단위는 라디안(Radian)입니다. 이러한 정보는 앞서 NCDF_LIST에 의하여 출력된 내용에서 변수 'x' 또는 'y'에 관한 부분에서 확인한 것입니다.

 

           4  x:                                                INT(5424) = INT(x)
                 0  scale_factor:   5.60000e-05
                 1  add_offset:     -0.151844
                 2  units: rad

 

다음으로 netCDF 파일로부터 읽어들일 다음 정보는 데이터 시각에 대한 변수 't'입니다. 먼저 NCDF_LIST에 의하여 출력된 관련 정보를 보면 다음과 같습니다.

 

           2  t:                                                DOUBLE  6.9549091e+08
                 0  long_name: J2000 epoch mid-point between the start and end image scan in seconds
                 1  standard_name: time
                 2  units: seconds since 2000-01-01 12:00:00

 

즉 파일 내에 수록된 데이터 시각 정보는 2000년 1월 1일 12시를 기점으로 한 초 단위의 값이라는 의미입니다. 따라서 변수 't'를 읽고 처리하는 과정은 다음과 같습니다.

 

NCDF_VARGET, id, 't', tsec
tj = JULDAY(1, 1, 2000, 12, 0, 0)+tsec/(24.*60*60)
tstr = STRING(tj, FORMAT='(C(CYI4.4, "/", CMoI2.2, "/", CDI2.2, " ", CHI2.2, ":", CMI2.2, ":", CSI2.2))')
PRINT, tstr

 

이와 같이 Julian Day 단위로 값들 얻은 후 우리가 알아볼 수 있도록 문자형 변수 tstr로 변환하였습니다. 출력된 내용은 다음과 같습니다.

 

2022/01/15 03:55:05

 

앞서 웹페이지에서 03시 50분 데이터 파일을 받은 것이긴 하지만 10분 단위이기 때문에 중간에 해당되는 03시 55분으로 시각 정보가 수록된 것으로 보입니다. 이제 마지막으로 netCDF 파일로부터 나머지 필요한 정보들을 다음과 같이 읽어들입니다.

 

NCDF_VARGET, id, 'nominal_satellite_subpoint_lon', clon
NCDF_VARGET, id, 'band_id', band_id
NCDF_VARGET, id, 'band_wavelength', band_wavelength
HELP, clon, band_id, band_wavelength

 

여기서는 중심 경도, 밴드 번호, 밴드 파장 등과 같은 정보들을 추출하였습니다. 출력된 정보는 다음과 같습니다.

 

CLON            FLOAT     =      -137.200
BAND_ID         BYTE      =   13
BAND_WAVELENGTH FLOAT     =       10.3300

 

이와 같이 추출된 변수 및 정보들을 이용하여 표출 과정까지 가보겠습니다. 이 과정에서는 GOES-R 투영법 기반의 지도를 표출하고 그 위에 이미지 데이터를 중첩하게 됩니다. 이 과정은 다음과 같습니다.

 

win = WINDOW(DIMENSIONS=[800, 700], /NO_TOOLBAR)
m = MAP('GOES-R', CENTER_LONGITUDE=clon, $
  HORIZON_COLOR='black', HORIZON_LINESTYLE=6, $
  CLIP=0, MARGIN=[0.05, 0.05, 0.15, 0.05], /CURRENT)
ct = COLORTABLE(57, /REVERSE)
im = IMAGE(img, x, y, MIN_VALUE=0, MAX_VALUE=130, $
  GRID_UNITS=1, RGB_TABLE=ct, /OVERPLOT)
mc = MAPCONTINENTS(/HIRES, COLOR='tomato')
m.MapGrid.LINESTYLE = 0
m.MapGrid.COLOR= 'black'
m.MapGrid.GRID_LONGITUDE = 15
m.MapGrid.GRID_LATITUDE = 15
cbar = COLORBAR(TARGET=im, ORIENTATION=1, POSITION=[0.94, 0.2, 0.98, 0.8], $
  TITLE='Radiance (mW m!e-2!n sr!e-1!n (cm!e-1!n)!e-1!n)', /BORDER)
tx0 = TEXT(0.01, 0.97, 'GOES-17 ABI Radiance (Band 13 : 10.35um)', FONT_SIZE=14, /NORMAL)
tx = TEXT(0.99, 0.97, tstr, ALIGNMENT=1, FONT_SIZE=14, /NORMAL)

 

여기서는 먼저 MAP 함수를 사용하여 GOES-R 투영법 기반의 바탕 지도를 표출한 다음 IMAGE 함수를 사용하여 이미지 데이터인 img를 중첩하는 것이 핵심입니다. 여기서 유의할 점은 이미지 배열인 img 뿐만 아니라 격자들의 X 및 Y 방향 위치값들에 해당되는 x 및 y 배열들도 함께 사용해야 한다는 것입니다. 그리고 이 위치값들은 scan-angle radian에 해당되는 거리값들이기 때문에 IMAGE 함수 내에서 GRID_UNITS 속성의 값이 1이 되어야 한다는 것도 유의해야 합니다. 제가 예전에 이 블로그를 통해서 천리안 2A호 위성 데이터를 GOES-R 투영법 지도상에 중첩 표출하는 방법에 관하여 소개했던 바 있는데, 기본적인 흐름은 이와 유사하다고 보시면 됩니다. 따라서 이 게시물의 내용도 함께 참조하시면 좋을 것 같습니다. 어쨌든 위의 과정에 의하여 표출된 결과는 다음 그림과 같습니다.

 

 

따라서 GOES 16/17 데이터를 netCDF 파일로 얻어서 이 파일을 읽고 데이터를 표출하는 과정은 대략 이러한 맥락이라고 보면 될 것 같습니다. 물론 데이터의 종류에 따라서 세부적인 차이들은 있을 것입니다.

 

그리고 오늘 소개한 내용에서 사용했던 예제 데이터의 날짜가 UTC로 2022년 1월 15일 03시 50분이었는데, 사실 이 때가 통가(Tonga) 화산 폭발이 있었던 시점인 04시 직전에 해당됩니다. 그래서 제가 03시 50분부터 11시 50분까지 7시간에 걸친 데이터들을 읽어서 동영상으로 만들어보았는데, 이 동영상은 아래 링크에서 보실 수 있습니다. 왼쪽 하단의 노란색 사각형이 바로 화산 폭발이 발생했던 부분입니다. 참고로 봐주시기 바랍니다.

 

https://youtu.be/Lz8F_juuzmA