NCL

NCL로 상층 등고선과 풍향풍속 그려보기

이석사 중 2023. 7. 20. 11:07
728x90

이번에는 850, 700, 500, 300hPa에서의 지위고도 선과 풍향 풍속을 그려보겠습니다


이번 코드에는 조건문이 들어갑니다

 

각 일기도마다 그려지는 자료들이 다르기 때문입니다

 

자주 쓴건 아니지만 어려운 부분은 아닙니다

 

주석으로 달아놓겠습니다

 

자료는 저번 포스팅에서 봤던 WRF 결과를 사용했습니다


코드입니다

begin
  f = addfile("./wrfout_d01_2016-10-06_03.nc", "r")
  times = wrf_user_list_times(f)
  ntimes = dimsizes(times)

  wks = gsn_open_wks("png", "wrf")

  pressure_levels = (/850., 700., 500., 300./) #압력 레벨

  nlevels = dimsizes(pressure_levels)

  hgt_spacing = (/20., 30., 60., 60./) #등고선 간격

  line_res1 = True
  line_res1@cnLineColor = "Red"
  line_res1@ContourParameters = (/5.0/)
  line_res1@cnInfoLabelOrthogonalPosF = 0.07
  line_res1@gsnContourLineThicknessesScale = 2.0

  fill_res = True #상대습도 resource
  fill_res@cnFillOn = True
  fill_res@ContourParameters = (/10., 90., 10./) 
  fill_res@cnFillColors = (/"White", "White", "White", \
  "White", "Chartreuse", "Green", \
  "Green3", "Green4", "ForestGreen", "PaleGreen4"/)

  line_res2 = True
  line_res2@cnLineColor = "MediumSeaGreen"
  line_res2@ContourParameters = (/10./)
  line_res2@cnInfoLabelOrthogonalPosF = 0.07
  line_res2@gsnContourLineThicknessesScale = 3.0

  vec_res = True
  vec_res@FieldTitle = "Wind"
  vec_res@NumVectors = 37

  line_res3 = True
  line_res3@cnLineColor = "NavyBlue"
  line_res3@gsnContourLineThicknessesScale = 3.0

  do it = 0, ntimes-1
    print("Working on time :" + times(it))
    line_res1@TimeLabel = times(it)
    line_res2@TimeLabel = times(it)
    fill_res@TimeLabel = times(it)

    tc = wrf_user_getvar(f, "tc", it)
    u = wrf_user_getvar(f, "ua", it)
    v = wrf_user_getvar(f, "va", it)
    p = wrf_user_getvar(f, "pressure", it)
    z = wrf_user_getvar(f, "z", it)
    rh = wrf_user_getvar(f, "rh", it)
    wspd_wdir = wrf_user_getvar(f, "wspd_wdir", it)

    do nl = 0, nlevels-1
      pressure = pressure_levels(nl)

      tc_plane = wrf_user_interp_level(tc, p, pressure, False)
      z_plane = wrf_user_interp_level(z, p, pressure, False)
      rh_plane = wrf_user_interp_level(rh, p, pressure, False)
      u_plane = wrf_user_interp_level(u, p, pressure, False)
      v_plane = wrf_user_interp_level(v, p, pressure, False)

      wspd_plane = wrf_user_interp_level(wspd_wdir(0, :, :, :), p, pressure, False)
      wdir_plane = wrf_user_interp_level(wspd_wdir(1, :, :, :), p, pressure, False)

      wspd_plane = (u_plane * u_plane + v_plane * v_plane)^(0.5)
      u_plane = u_plane * 1.94386
      v_plane = v_plane * 1.94386
      wspd_plane@description = "Wind Speed"
      wspd_plane@units = "m/s"
      u_plane@units = "kts"
      v_plane@units = "kts"

      line_res3@ContourParameters = hgt_spacing(nl)

      contour_tc = wrf_contour(f, wks, tc_plane, line_res1)
      contour_rh = wrf_contour(f, wks, rh_plane, fill_res)
      contour_wspd = wrf_contour(f, wks, wspd_plane, line_res2)
      vector = wrf_vector(f, wks, u_plane, v_plane, vec_res)
      contour_height = wrf_contour(f, wks, z_plane, line_res3)

;고도별로 입력되는 자료들이 다르기 때문에 조건문으로 고도별로 그림
      if (pressure .eq. 850) then 
        ov = wrf_map_overlays(f, wks, (/contour_rh, contour_tc, contour_height, vector/), True, True)
      elseif (pressure .eq. 700) then
        ov = wrf_map_overlays(f, wks, (/contour_tc, contour_height, vector/), True, True)
      elseif (pressure .eq. 500) then
        ov = wrf_map_overlays(f, wks, (/contour_tc, contour_height, vector/), True, True)
      elseif (pressure .eq. 300) then
        ov = wrf_map_overlays(f, wks, (/contour_wspd, contour_height, vector/), True, True)
      end if
    end do
  end do
end

결과 입니다

 

결과는 총 4장이 그려집니다

728x90