[IDL] 아이디엘 대기 환경 자료 처리 및 가시화 : Grib 형식인 ECMWF Interim 기상 자료를 이용한 가시화

정보
내용
소스 코드
참고 문헌
문의사항

정보

업무명 : 대기 환경 자료 처리 및 가시화 : Grib 형식인 ECMWF Interim기상 자료를 이용한 가시화
작성자 : 이상호
작성일 : 2020-03-31
설 명 :
수정이력 :

내용

[개요]

안녕하세요? 기상 연구 및 웹 개발을 담당하고 있는 해솔입니다.
대기 환경 정보는 다양한 차원의 관측 자료 및 분석 자료로 구성되어 있기 때문에 이용자가 적절한 분석 및 가시화 기술이 요구됩니다.
따라서 오늘 포스팅에서는 Grib 형식인 ECMWF Interim 기상 자료를 이용한 가시화를 소개해 드리고자 합니다.
추가로 대기과학 전공자를 위한 IDL를 소개한 링크를 보내드립니다.

[IDL] 아이디엘 대기과학 전공자를 위한 IDL (Interactive Data Language) 소개

정보 업무명 : 대기과학 전공자를 위한 IDL (Interactive Data Language) 소개 작성자 : 이상호 작성일 : 2020-03-23 설 명 : 수정이력 : 내용 [개요] 안녕하세요? 웹 개발 및 연구 개발을 담당하고 있는 해솔입니..

shlee1990.tistory.com

[특징]

대기 환경 자료를 이해하기 위해서 자료 처리 및 가시화 기술이 요구되며 이 프로그램은 이러한 목적을 달성하기 위해 고안된 소프트웨어

[기능]

위도, 경도, ECMWF Interim 자료 읽기
자료 전처리
가시화

[활용 자료]

자료명 : interim_sfc_2014102500_00.grib
자료 종류 : 기상 자료 (지표면 기압, 지표면 온도 등)
확장자 : grib
영역 : 전지구
기간 : 2014년 10월 25일

interim_sfc_2014102500_00.z01

10.00MB

interim_sfc_2014102500_00.zip

5.97MB

[사용법]

입력 자료를 동일 디렉터리에 위치
소스 코드를 실행 (idl -e Visualization_Using_Weather_Data_of_Grib_Format)
가시화 결과를 확인

[사용 OS]

Windows 10

[사용 언어]

IDL v8.5

소스 코드

해당 작업을 하기 위한 컴파일 및 실행 코드는 IDL로 작성되었으며 가시화를 위한 라이브러리는 Coyote's Guide to IDL Programming를 이용하였습니다.
소스 코드는 단계별로 수행하며 상세 설명은 다음과 같습니다.
- 1 단계는 주 프로그램은 작업 경로 설정, Grib 파일 읽기, 변수 전처리하여 메모리상에 저장하고 가시화를 위한 초기 설정합니다.
- 2 단계는 plot 매핑에 따라 영상 장면 표출하여 이미지 형식으로 저장합니다. 이 과정에서 포스트 스크립트 (PS) 형식에서 PNG로 변환합니다.

[명세]

전역 변수 설정
- cd를 통해 작업 디렉터리 설정
- grib 파일명 설정 (gribFile)
- 디코딩을 위한 변수명 설정 (varName)
- 배열 정보 정의 (dim)
- 위도 (lat1D), 경도 (lon1D), 자료값 (val2D)에 대한 배열 설정

cd, 'C:\SYSTEM\PROG\IDL\GRIB'

gribFile = 'INPUT/interim_sfc_2014102500_00.grib'

varName = ['srp', 'tcc', '10u', '10v', '2mT', 'lcc', 'mcc', 'hcc', 'skt']

dim = [1440, 721]

lon1D = fltarr(dim[0])
lat1D = fltarr(dim[1])
val2D = fltarr(dim[0], dim[1])

위/경도 설정
- 위도 및 경도는 각각 25 km 간격으로 1440 및 721개 격자 구성
- 경도가 180 도 이상일 경우 -360 도를 추가

lon1D = 0.0 + FINDGEN(dim[0]) * 0.25
lat1D = 90.0 - FINDGEN(dim[1]) * 0.25

loc1 = where(lon1D gt 180.0, nloc1)
if (nloc1 gt 0) then lon1D[loc1] = lon1D[loc1] - 360

Grib 파일에서 변수 처리 및 파일 읽기
- wgrib.exe를 통해 변수 (varName)마다 출력 경로 (outFile)에 저장
- 출력 파일 (outFile)을 읽기
- 추가로 변수 디코딩에 필요한 wgrib.exe 설치 방법에 대한 링크를 보내드립니다.

[Python] 파이썬 원도우 (Window 10)에서 수치예측모델 (Grib, Grb2) 자료 처리를 위한 "wgrib, wgrib2, cdo" 설치 방법

정보 업무명 : 파이썬 원도우 (Window 10)에서 수치예측모델 (Grib, Grb2) 자료 처리를 위한 wgrib, wgrib2, cdo 설치 방법 작성자 : 이상호 작성일 : 2020-02-09 설 명 : 수정이력 : 내용 [개요] 안녕하세요? 웹..

shlee1990.tistory.com

for iCount = 0L, N_ELEMENTS(varName) - 1 do begin

    varId = strtrim(string(iCount + 1), 2)

    outFile = './OUTPUT/' + varName(iCount) + '.txt'

    com = 'wgrib.exe ' + gribFile + ' -i -d ' + varId + ' -text -nh -o ' + outFile     

    spawn, /hide, com

    close, 1
    openr, 1, outFile & readf, 1, val2D & close, 1

	# 중략
    
endfor

요약 통계량 확인
- 최소값 (min), 평균값 (mean), 최대값 (max) 계산
- help를 통해 자료형 및 배열 정보 확인

print, min(lon1D, /nan), mean(lon1D, /nan), max(lon1D, /nan)
print, min(lat1D, /nan), mean(lat1D, /nan), max(lat1D, /nan)
print, min(val2D, /nan), mean(val2D, /nan), max(val2D, /nan)

help, lon1D, lat1D, val2D

가시화
- 데이터 세트에서 정보는 여전히 숨겨져 있기 때문에 가시화 필요
- 가시화는 일반적인 정적 탐색 데이터 분석에서 웹 브라우저의 동적 대화식 데이터 시각화에 이르기까지 다양함

Grib 변수를 위한 설정 및 가시화
- 관심 영역 (ROI)와 중심 위/경도 (cnLon, cnLat) 설정
- 이미지 저장 파일명 (psName) 및 제목명 (mainName)
- switch 및 fnMakePsPlot를 통해 지표면 기압과 지표면 온도 이미지 저장

ROI = [-90, 90, -180, 180]

cnLat = 0
cnLon = 0

dtDateYmdH = STRMID(Gribfile, 18, 13)

switch iCount of

    0 : begin 
        psName = "Img01_" + dtDateYmdH
        mainName = "Grib Surface Pressure : " + dtDateYmdH

        fnMakePsPlot, lon1D, lat1D, val2D / 100.0, dim, cnLat, cnLon, ROI, psName, mainName, 600, 1200, 200, 33
    end

    4 : begin 
        psName = "Img02_" + dtDateYmdH
        mainName = "Grib Surface Air Temperature : " + dtDateYmdH

        fnMakePsPlot, lon1D, lat1D, val2D, dim, cnLat, cnLon, ROI, psName, mainName, 250, 320, 10, 33
    end
endswitch

지표면 기압

지표면 온도

사용자 편의성 함수 정의
- 전달 인자
  - 1차원 위도 (lon1D), 1차원 경도 (lat1D), 2차원 값 (val2D), 배열 정보 (dim), 중심 위도 (cnLat), 중심 경도 (cnLon), 관심 영역 (ROI), 이미지 저장 파일명 (psName), 제목명 (mainName), 최소값 (zmin), 최대값 (zmax), 등고선 간격 (interval), 컬러 테이블 (color_table)
- 기능
  - plot를 통해 위/경도에 따른 값 매핑
  - 포스트 스크립트 (PS)를 PNG로 변환

;====================================================================
; Subroutine : Function Make Postscript to Png Graph
;====================================================================
pro fnMakePsPlot, lon1D, lat1D, val2D, dim, cnLat, cnLon, ROI, psName, mainName, zmin, zmax, interval, color_table

    set_plot, "ps"

    device, filename = psName + ".ps", decomposed=0, bits=8, /color, xsize=16, ysize = 12, /inches, font_size = 13, /Helvetica

    !p.font = 0 & !p.charsize = 2.0 & !p.charthick = 1.6  & !p.multi = [0,1,1] & !p.background = 255

    xvert = [0.1, 0.1, -0.1, -0.1 ,0.1]
    yvert = [-0.1, 0.1, 0.1, -0.1, -0.1]

    usersym, xvert, yvert, /fill

    start_color = 0 & end_color = 255 & colorn = end_color - start_color + 1

    cgloadct, color_table

    latmin = ROI(0) & latmax= ROI(1) & lonmin = ROI(2) & lonmax= ROI(3)

    cgMAP_SET, cnLat, cnLon, /CYLINDRICAL $
      , limit = [latmin, lonmin, latmax, lonmax] $
      , position = [0.1, 0.1, 0.8, 0.85], /noborder

    for i = 0L, dim[0] - 1 do begin
      for j = 0L, dim[1] - 1 do begin

        nLon = lon1D[i]
        nLat = lat1D[j]
        nVal = val2D[i, j]

        if (lonmin gt nLon or nLon gt lonmax or finite(nLon, /nan) eq 1) then continue
        if (latmin gt nLat or nLat gt latmax or finite(nLat, /nan) eq 1) then continue
        if (nVal le 0 or finite(nVal, /nan) eq 1) then continue

        plots, nLon, nLat, psym = 8, symsize = 2, color = BYTSCL(nVal, zmin, zmax)

      endfor
    endfor

    cgColorbar, NColors = 255, BOTTOM = 1 $
      , Divisions = 5, COLOR = "black" $
      , Position = [0.2, 0.90, 0.8, 0.94] $
      , Range = [zmin, zmax] $
      , CHARSIZE = 1.5,CHARthick = 2 $
      , VERTICAL = 1, RIGHT = 1

      cgLoadct, 2

      contour, val2D, lon1D, lat1D, /overplot, C_THICK = 2 $
        , LEVELS = zmin + findgen(10) * interval $
        , C_LABELS = zmin + findgen(10) * interval $
        , C_CHARSIZE = 1.5, C_CHARTHICK = 2 $
        , C_COLORS = BYTSCL(zmin + findgen(10) * interval, zmin, zmax)

    cgMAP_CONTINENTS, /coast, /countries, COLOR = "black"

    lats = [-90:90:30]
    lons = [-180:360:60]

    lats_names = strarr(n_elements(lats))
    lons_names = strarr(n_elements(lons))

    for i = 0L, n_elements(lats) - 1 do begin
      if (lats[i] gt 0) then begin
        lats_names[i] = textoidl(string(lats[i]) + "\circN")
      endif else if (lats[i] eq 0) then begin
        lats_names[i] = textoidl(string(lats[i]) + "\circ")
      endif else begin
        lats_names[i] = textoidl(string(-lats[i]) + "\circS")
      endelse
    endfor

    for i = 0L,n_elements(lons) - 1  do begin
      if (lons[i] gt 0) then begin
        lons_names[i] = textoidl(string(lons[i]) + "\circE")
      endif else if (lons[i] eq 0) then begin
        lons_names[i] = textoidl(string(lons[i]) + "\circ")
      endif else begin
        lons_names[i] = textoidl(string(-lons[i]) + "\circW")
      endelse
    endfor

    cgMap_grid, color = "black", charsize = 1.8, thick = 2.5, lats = lats,  latnames = lats_names, lons = lons, lonnames = lons_names $
      , linestyle = 1, bthick = 300, /box_axes, /no_grid

    cgText, 0.45, 0.95, mainName, /normal, charsize = 2, CHARTHICK = 2, color = "black", alignment = 0.5

    device, /close_file

    com = "convert -flatten -background white " + psName + ".ps" + " " + "./FIG/" + file_basename(psName, ".ps") + ".png"
    spawn, /hide, com

    file_delete, psName + ".ps"

    return

end

[전체]

	pro Visualization_Using_Weather_Data_of_Grib_Format

	cd, 'C:\SYSTEM\PROG\IDL\GRIB'

	gribFile = 'INPUT/interim_sfc_2014102500_00.grib'

	varName = ['srp', 'tcc', '10u', '10v', '2mT', 'lcc', 'mcc', 'hcc', 'skt']

	dim = [1440, 721]

	lon1D = fltarr(dim[0])
	lat1D = fltarr(dim[1])
	val2D = fltarr(dim[0], dim[1])

	lon1D = 0.0 + FINDGEN(dim[0]) * 0.25
	lat1D = 90.0 - FINDGEN(dim[1]) * 0.25

	loc1 = where(lon1D gt 180.0, nloc1)
	if (nloc1 gt 0) then lon1D[loc1] = lon1D[loc1] - 360


	for iCount = 0L, N_ELEMENTS(varName) - 1 do begin

	varId = strtrim(string(iCount + 1), 2)

	outFile = './OUTPUT/' + varName(iCount) + '.txt'

	; print, varId, varName(iCount)

	com = 'wgrib.exe ' + gribFile + ' -i -d ' + varId + ' -text -nh -o ' + outFile

	spawn, /hide, com

	close, 1
	openr, 1, outFile & readf, 1, val2D & close, 1


	print, min(lon1D, /nan), mean(lon1D, /nan), max(lon1D, /nan)
	print, min(lat1D, /nan), mean(lat1D, /nan), max(lat1D, /nan)
	print, min(val2D, /nan), mean(val2D, /nan), max(val2D, /nan)

	help, lon1D, lat1D, val2D


	ROI = [-90, 90, -180, 180]

	cnLat = 0
	cnLon = 0

	dtDateYmdH = STRMID(Gribfile, 18, 13)

	switch iCount of

	0 : begin
	psName = "Img01_" + dtDateYmdH
	mainName = "Grib Surface Pressure : " + dtDateYmdH

	fnMakePsPlot, lon1D, lat1D, val2D / 100.0, dim, cnLat, cnLon, ROI, psName, mainName, 600, 1200, 200, 33
	end

	4 : begin
	psName = "Img02_" + dtDateYmdH
	mainName = "Grib Surface Air Temperature : " + dtDateYmdH

	fnMakePsPlot, lon1D, lat1D, val2D, dim, cnLat, cnLon, ROI, psName, mainName, 250, 320, 10, 33
	end
	endswitch

	endfor

	end

	;====================================================================
	; Subroutine : Function Make Postscript to Png Graph
	;====================================================================
	pro fnMakePsPlot, lon1D, lat1D, val2D, dim, cnLat, cnLon, ROI, psName, mainName, zmin, zmax, interval, color_table

	set_plot, "ps"

	device, filename = psName + ".ps", decomposed=0, bits=8, /color, xsize=16, ysize = 12, /inches, font_size = 13, /Helvetica

	!p.font = 0 & !p.charsize = 2.0 & !p.charthick = 1.6 & !p.multi = [0,1,1] & !p.background = 255

	xvert = [0.1, 0.1, -0.1, -0.1 ,0.1]
	yvert = [-0.1, 0.1, 0.1, -0.1, -0.1]

	usersym, xvert, yvert, /fill

	start_color = 0 & end_color = 255 & colorn = end_color - start_color + 1

	cgloadct, color_table

	latmin = ROI(0) & latmax= ROI(1) & lonmin = ROI(2) & lonmax= ROI(3)

	cgMAP_SET, cnLat, cnLon, /CYLINDRICAL $
	, limit = [latmin, lonmin, latmax, lonmax] $
	, position = [0.1, 0.1, 0.8, 0.85], /noborder

	for i = 0L, dim[0] - 1 do begin
	for j = 0L, dim[1] - 1 do begin

	nLon = lon1D[i]
	nLat = lat1D[j]
	nVal = val2D[i, j]

	if (lonmin gt nLon or nLon gt lonmax or finite(nLon, /nan) eq 1) then continue
	if (latmin gt nLat or nLat gt latmax or finite(nLat, /nan) eq 1) then continue
	if (nVal le 0 or finite(nVal, /nan) eq 1) then continue

	plots, nLon, nLat, psym = 8, symsize = 2, color = BYTSCL(nVal, zmin, zmax)

	endfor
	endfor

	cgColorbar, NColors = 255, BOTTOM = 1 $
	, Divisions = 5, COLOR = "black" $
	, Position = [0.2, 0.90, 0.8, 0.94] $
	, Range = [zmin, zmax] $
	, CHARSIZE = 1.5,CHARthick = 2 $
	, VERTICAL = 1, RIGHT = 1

	cgLoadct, 2

	contour, val2D, lon1D, lat1D, /overplot, C_THICK = 2 $
	, LEVELS = zmin + findgen(10) * interval $
	, C_LABELS = zmin + findgen(10) * interval $
	, C_CHARSIZE = 1.5, C_CHARTHICK = 2 $
	, C_COLORS = BYTSCL(zmin + findgen(10) * interval, zmin, zmax)

	cgMAP_CONTINENTS, /coast, /countries, COLOR = "black"

	lats = [-90:90:30]
	lons = [-180:360:60]

	lats_names = strarr(n_elements(lats))
	lons_names = strarr(n_elements(lons))

	for i = 0L, n_elements(lats) - 1 do begin
	if (lats[i] gt 0) then begin
	lats_names[i] = textoidl(string(lats[i]) + "\circN")
	endif else if (lats[i] eq 0) then begin
	lats_names[i] = textoidl(string(lats[i]) + "\circ")
	endif else begin
	lats_names[i] = textoidl(string(-lats[i]) + "\circS")
	endelse
	endfor

	for i = 0L,n_elements(lons) - 1 do begin
	if (lons[i] gt 0) then begin
	lons_names[i] = textoidl(string(lons[i]) + "\circE")
	endif else if (lons[i] eq 0) then begin
	lons_names[i] = textoidl(string(lons[i]) + "\circ")
	endif else begin
	lons_names[i] = textoidl(string(-lons[i]) + "\circW")
	endelse
	endfor

	cgMap_grid, color = "black", charsize = 1.8, thick = 2.5, lats = lats, latnames = lats_names, lons = lons, lonnames = lons_names $
	, linestyle = 1, bthick = 300, /box_axes, /no_grid

	cgText, 0.45, 0.95, mainName, /normal, charsize = 2, CHARTHICK = 2, color = "black", alignment = 0.5

	device, /close_file

	com = "convert -flatten -background white " + psName + ".ps" + " " + "./FIG/" + file_basename(psName, ".ps") + ".png"
	spawn, /hide, com

	file_delete, psName + ".ps"

	return

	end

view raw Visualization_Using_Weather_Data_of_Grib_Format.pro hosted with ❤ by GitHub

참고 문헌

[논문]

없음

[보고서]

없음

[URL]

없음

문의사항

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

sangho.lee.1990@gmail.com

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

saimang0804@gmail.com

저작자표시 비영리 변경금지

'프로그래밍 언어 > IDL' 카테고리의 다른 글

[ENVI] ENVI를 이용한 원격탐사 (Remote Sensing with ENVI) 교육 연수 (0)	2020.04.03
[IDL] 아이디엘 대기 환경 자료 처리 및 가시화 : 태양 위치 정보를 이용하여 태양 천정각 (Solar Zenith Angle) 및 태양 방위각 (Solar Azimuth Angle) 가시화 (0)	2020.04.01
[IDL] 아이디엘 대기 환경 자료 처리 및 가시화 : HDF5 형식인 천리안위성 1A호 (COMS/MI) 기상 위성 자료를 이용한 에어로졸 광학두께 가시화 (4)	2020.03.31
[IDL] 아이디엘 대기 환경 자료 처리 및 가시화 : HDF 형식인 전지구 오존량 자료를 이용한 가시화 (0)	2020.03.31
[IDL] 아이디엘 대기 환경 자료 처리 및 가시화 : NetCDF 형식인 대기 온도 자료를 이용한 가시화 (0)	2020.03.31

내 블로그 - 관리자 홈 전환	`Q` `Q`
새 글 쓰기	`W` `W`

글 수정 (권한 있는 경우)	`E` `E`
댓글 영역으로 이동	`C` `C`

이 페이지의 URL 복사	`S` `S`
맨 위로 이동	`T` `T`
티스토리 홈 이동	`H` `H`
단축키 안내	`Shift` + `/` `⇧` + `/`