[Python] 파이썬 바이너리 (Binary) 형식인 천리안위성 1A호 (COMS/MI) 기상 위성 자료를 이용한 가시화

 정보

• 업무명     : 바이너리 (Binary) 형태인 천리안위성 1A호 (COMS/MI) 기상 위성 자료를 이용한 가시화

• 작성자     : 이상호

• 작성일     : 2019-09-03

• 설   명      :

• 수정이력 :

  • 2020-02-22 : 소스 코드 명세 추가

 

 내용

 [개요]

• 안녕하세요? 웹 개발 및 연구 개발을 담당하고 있는 해솔입니다.

• 대기과학에서 제공되는 자료 형식은 "CSV, 바이너리 이진 자료 (Binary), NetCDF, HDF4, HDF5, Grib, Grb2, Bufr" 등이 있습니다.

• 대부분 용도에 맞게 사용하기 때문에 각각 모든 형식에 맞게 자료 전처리 및 가시화는 필연입니다.

• 예를 들면 기상청의 경우 지상 관측 자료와 기상 위성 자료 및 수치예보 그리고 기타의 경우 각각 "CSV 형식"와 "NetCDF, HDF4, HDF5" 및 "Grib, Gr2, Bufr" 그리고 "바이너지 이진 자료 (Binary)"로 배포합니다.

• 오늘은 대기과학에서 사용되는 바이너리 이진 자료 (Binary) 설명뿐만 아니라 Python을 이용하여 자료 전처리 및 가시화를 소개해 드리고자 합니다. 

 

 

[특징]

• 바이너리 형태인 기상위성 자료를 이해하기 위해 가시화 도구가 필요하며 이 프로그램은 이러한 목적을 달성하기 위해 고안된 소프트웨어

 

[기능]

• 천리안위성 1호 (COMS/MI) 기상위성 자료를 통해 도법에 따른 가시화

 

[활용 자료]

• 위성명 : 천리안위성 1A호 (COMS) 기상위성 자료

• 센서명 : MI 기상영상기

• 자료종류 : 대기 상단에서의 반사율

• 영역 : 전구

• 해상도 : 5 km

• 도법 : Lambert Conformal Conic 및 Orthographic

• 확장자 : 바이너리 (.bin)

• 기간 : 2017년 08월 04일 0000 UTC (0900 KST)

 

[자료 처리 방안 및 활용 분석 기법]

• 없음

 

[사용법]

• 입력 자료를 동일 디렉터리 위치

• 소스 코드를 실행 (Jupyter notebook 실행)

• 가시화 결과를 확인

 

[사용 OS]

• Window 10

 

[사용 언어]

• Python v2.7

 

 소스 코드

[명세]

• 라이브러리 읽기

# Library
import pandas as pd
import numpy as np
import sys 
import os
import matplotlib.pyplot as plt
from dplython import *
    # (DplyFrame, X, diamonds, select, sift, sample_n,
    #     sample_frac, head, arrange, mutate, group_by, summarize, DelayFunction) 
from scipy.stats import linregress
from matplotlib import pyplot as plt
from IPython.display import Image
from mpl_toolkits.basemap import Basemap
from matplotlib.colors import Normalize

import matplotlib
import matplotlib.cm as cm
import seaborn as sns
from scipy.stats import linregress
from matplotlib import rcParams
#from netCDF4 import Dataset # NC 파일을 읽기 위한 사전작업 실시
import struct
import binascii
from pyhdf.SD import SD, SDC # HDF4 읽는 라이브러리
import h5py
import cython
#import ccplot.utils
#from ccplot.hdf import *

 

• 위/경도 이진자료 읽기

  • 격자수 (1950000 * 4) 및 32 bit 실수형 (f4)

  • 1차원 위/경도 변수

nx, ny = 1950000, 4
latlon = np.fromfile('MI/LatLon/East_Asia.bin', dtype='f4', count=nx*ny)

print(latlon.shape)

 

• 위/경도에 대해 추출

latlon_L1 = np.reshape(latlon, (-1, 4))

lat1D = latlon_L1[:,2]
lon1D = latlon_L1[:,3]

 

• 가시 채널 이진자료 읽기

  • 격자수 (1950000) 및 자료형 (big-endian unsigned short)

nx = 1950000
dn1D = np.fromfile('MI/coms_le1b_vis_ch5_a_201708040000.bin', \
                  dtype='>H', count=nx)
print(dn1D.shape)

 

• Data Frame 설정

  • 1차원 위도, 경도, 변수 (rsr1D, dsr1D, dn1D)

data =  pd.DataFrame( np.column_stack( [lat1D, lon1D, dn1D] ),
                     columns=['lat', 'lon', 'dn'], )
data.head()

 

	lat	lon	dn
0	46.726562	78.285103	69.0
1	46.756237	78.323036	69.0
2	46.785904	78.361023	72.0
3	46.815556	78.399071	71.0
4	46.845200	78.437172	72.0

 

• Data Frame를 통해 L1 전처리

  • 최근 데이터 분석에서 사용되는 R의 "dplyr" 라이브러리와 유사한 "dplython"를 사용

  • DN 변수에 대한 최대값 및 최소값 설정

  • "COMS/MI에 대한 조견표 (Look-Up-Table)"를 참조하여 DN 값을 복사량 (Radiance) 및 알벤도 (albedo) 또는 밝기 온도 (Brightness Temperature)

data_L1 = ( DplyFrame(data) >>
           sift( (50 <= X.dn) & (X.dn <= 901) ) >>
           mutate( radiance = (0.6463*X.dn) - 32.0544 ) >>
           mutate( albedo   = (0.2091*X.dn) + 1.7014E-005 )        
          )

data_L1.head()

 

lat	lon	dn	radiance	albedo
0	46.726562	78.285103	69.0	12.5403	14.427917
1	46.756237	78.323036	69.0	12.5403	14.427917
2	46.785904	78.361023	72.0	14.4792	15.055217
3	46.815556	78.399071	71.0	13.8329	14.846117
4	46.845200	78.437172	72.0	14.4792	15.055217

 

• 가시화를 위한 설정

• 초기 설정

  • 크기 : figure

  • 스타일 : style

  • 폰트 : rc, rcParams

  • 지도 해상도 : Basemap

     

• 변수 설정

  • 위도 (lon), 경도 (lat), 알베도 (albedo)

 

• 그림 설정

  • 산점도 : scatter

  • 컬러바 : colorbar

  • 해안선 : drawcostlines, drawcountries, drawmapboundary

  • 수평/수직 그리드 : drawparallels, drawmeridians

  • 그림 제목 : title

  • 그림 위치에 따른 텍스트 삽입 : annotate

 

• Lambert Conformal Conic 도법

%matplotlib inline
# define plot size in inches (width, height) & resolution(DPI)
plt.figure(figsize=(12, 10))

# style
plt.style.use('seaborn-darkgrid')

# define font size
plt.rc("font", size=22)
plt.rcParams['font.family'] = 'New Century Schoolbook'
# plt.rcParams["font.weight"] = "bold"
# plt.rcParams['axes.labelsize'] = 26
# plt.rcParams['xtick.labelsize'] = 26
# plt.rcParams['ytick.labelsize'] = 26
# plt.rcParams["axes.labelweight"] = "bold"
# plt.rcParams["axes.titleweight"] = "bold"

m = Basemap(projection='lcc', lon_0=128.2, lat_0=0, lat_1=10, lat_2=10, 
            llcrnrlon=75, llcrnrlat=0, 
            urcrnrlon=190, urcrnrlat=60, 
            resolution='c')

# m = Basemap(projection='lcc', lon_0=128.2, lat_0=0, lat_1=25.,lat_2=40., 
#             llcrnrlon=100, urcrnrlon=170, llcrnrlat=5, urcrnrlat=55, resolution='c')
# m = Basemap(projection='lcc', lon_0=140.714793, lat_0=-0.007606, llcrnrlon=min(data_q.lon)-10, urcrnrlon=max(data_q.lon)+10, llcrnrlat=10, urcrnrlat=max(data_q.lat)+5, resolution='l')
# m = Basemap(projection='mill', lon_0=140.714793, lat_0=-0.007606, resolution='c')
# c (crude) < l (low) < i (intermediate) < h (high) < f (full)
# Source: http://seesaawiki.jp/met-python/d/Basemap

X, Y = m(data_L1.lon.values, data_L1.lat.values)
VAL = data_L1.albedo.values

# m.scatter(X, Y, c=VAL, s=1.0, marker="s", zorder=1, vmin=0, vmax=100, cmap=plt.cm.get_cmap('gist_gray'), alpha=1.0) 
m.scatter(X, Y, c=VAL, s=1.0, marker="s", zorder=1, vmin=0, vmax=100, cmap=plt.cm.get_cmap('rainbow'), alpha=1.0) 
# m.colorbar(location='right', label='Reflectance  [%]')
m.colorbar(location='bottom', label='Reflectance  [%]', pad=0.5)

# VAL = data_L1.radiance.values
# m.scatter(X, Y, c=VAL, s=1.0, marker="s", zorder=1, vmin=0, vmax=300, cmap=plt.cm.get_cmap('gist_gray'), alpha=1.0) 
# m.colorbar(location='right', label='Radiance  $\mathregular{[Wm^{-2}]}$')
# m.colorbar(location='bottom', label='Radiance  $\mathregular{[Wm^{-2}]}$', pad=0.5)

# Black background
# m.drawcoastlines(color='yellow')
# m.drawcountries(color='yellow')
# m.drawmapboundary(fill_color='black')
# m.drawparallels(np.arange(-150, 120, 15), labels=[1,0,0,0], dashes=[2,2], color='yellow')
# m.drawmeridians(np.arange(-180, 180, 20), labels=[0,0,0,1], dashes=[2,2], color='yellow')

# White background
m.drawcoastlines(color='black')
m.drawcountries(color='black')
m.drawmapboundary(fill_color='white')
m.drawparallels(np.arange(-15, 120, 15), labels=[1,0,0,0], dashes=[2,2], color='black')
m.drawmeridians(np.arange(-180, 180, 15), labels=[0,0,0,1], dashes=[2,2], color='black')

plt.title('COMS / MI  2017-08-04  0000 UTC \n')
plt.show()

 

그림. 천리안위성 1A호 (COMS/MI) 기상위성 자료를 이용한 대기상단에서의 반사율 맵핑.

 

• Orthographic 도법

%matplotlib inline
# define plot size in inches (width, height) & resolution(DPI)
plt.figure(figsize=(12, 10))

# style
plt.style.use('seaborn-darkgrid')

# define font size
plt.rc("font", size=22)
plt.rcParams['font.family'] = 'New Century Schoolbook'
# plt.rcParams["font.family"] = "Century Schoolbook"
# rcParams['font.family'] = 'sans-serif'
# plt.rcParams["font.weight"] = "bold"
# plt.rcParams['axes.labelsize'] = 26
# plt.rcParams['xtick.labelsize'] = 26
# plt.rcParams['ytick.labelsize'] = 26
# plt.rcParams["axes.labelweight"] = "bold"
# plt.rcParams["axes.titleweight"] = "bold"

# m = Basemap(projection='ortho', lon_0=128.2, lat_0=0, resolution='c')
m = Basemap(projection='ortho', lon_0=140.714793, lat_0=-0.007606, resolution='c')

X, Y = m(data_L1.lon.values, data_L1.lat.values)
VAL = data_L1.albedo.values

m.scatter(X, Y, c=VAL, s=1.0, marker="s", zorder=1, vmin=0, vmax=100, cmap=plt.cm.get_cmap('rainbow'), alpha=1.0) 
# m.scatter(X, Y, c=VAL, s=1.0, marker="s", zorder=1, vmin=0, vmax=100, cmap=plt.cm.get_cmap('gist_gray'), alpha=1.0) 
m.colorbar(location='right', label='Reflectance  [%]')
# m.colorbar(location='bottom', label='Reflectance  [%]', pad=0.5)

m.drawcoastlines()
m.drawcountries()
m.drawmapboundary(fill_color='white')
m.drawparallels(np.arange(-60.,61.,30.),labels=[1,0,0,0],dashes=[2,2])
m.drawmeridians(np.arange(-160.,200.,30.),labels=[0,0,0,1],dashes=[2,2])
plt.title('COMS / MI  2017-08-04  0000 UTC \n')

plt.annotate(u'60\N{DEGREE SIGN}S', xy=(m(170, -62.5)), color='black', fontweight='bold', xycoords='data', horizontalalignment='center', verticalalignment='top')
plt.annotate(u'30\N{DEGREE SIGN}S', xy=(m(170, -32.5)), color='black', fontweight='bold', xycoords='data', horizontalalignment='center', verticalalignment='top')
plt.annotate(u'30\N{DEGREE SIGN}N', xy=(m(170, +27.5)), color='black', fontweight='bold', xycoords='data', horizontalalignment='center', verticalalignment='top')
plt.annotate(u'60\N{DEGREE SIGN}N', xy=(m(170, +57.5)), color='black', fontweight='bold', xycoords='data', horizontalalignment='center', verticalalignment='top')
plt.annotate(u'80\N{DEGREE SIGN}E', xy=(m(80, -2.5))  , color='black', fontweight='bold', xycoords='data', horizontalalignment='center', verticalalignment='top')
plt.annotate(u'110\N{DEGREE SIGN}E', xy=(m(110, -2.5)), color='black', fontweight='bold', xycoords='data', horizontalalignment='center', verticalalignment='top')
plt.annotate(u'110\N{DEGREE SIGN}E', xy=(m(110, -2.5)), color='black', fontweight='bold', xycoords='data', horizontalalignment='center', verticalalignment='top')
plt.annotate(u'140\N{DEGREE SIGN}E', xy=(m(140, -2.5)), color='black', fontweight='bold', xycoords='data', horizontalalignment='center', verticalalignment='top')
plt.annotate(u'170\N{DEGREE SIGN}E', xy=(m(170, -2.5)), color='black', fontweight='bold', xycoords='data', horizontalalignment='center', verticalalignment='top')
plt.annotate(u'160\N{DEGREE SIGN}W', xy=(m(-160, -2.5)),color='black', fontweight='bold', xycoords='data', horizontalalignment='center', verticalalignment='top')
# plt.savefig("FIG/Scane_analysis_Aqua_CERES_RSR.png")
plt.show()

 

그림. 천리안위성 1A호 (COMS/MI) 기상위성 자료를 이용한 대기상단에서의 반사율 가시화.

 

[전체]

 

 참고 문헌

[논문]

• 없음

[보고서]

• 없음

[URL]

• 없음

 

 문의사항

[기상학/프로그래밍 언어]

• sangho.lee.1990@gmail.com

[해양학/천문학/빅데이터]

• saimang0804@gmail.com