#----------------------------------------------------------------------
# Read and plots GRIB2 data in its native projection using PyNIO/PyNGL.
#
# Curvilinear grid, 461 x 421
#----------------------------------------------------------------------
import numpy as np
import Nio, Ngl

#---Read data
dir      = "../Data/"
filename = "MET9_IR108_cosmode_0909210000.grb2"
a    = Nio.open_file(dir+filename)
sbt  = a.variables["SBTMP_P31_GRLL0"][:]
lat  = a.variables["gridlat_0"]
lon  = a.variables["gridlon_0"]
nlat = sbt.shape[0]
nlon = sbt.shape[1]

wks = Ngl.open_wks("png","grib2_plot_rotated")

#---Set some plot options
res                   = Ngl.Resources()

# Contour options
res.cnFillOn           = True          # turn on contour fill
res.cnLinesOn          = False         # turn off contour lines
res.cnLineLabelsOn     = False         # turn off line labels
res.cnFillPalette      = "BlGrYeOrReVi200"
res.cnLevelSpacingF    = 10

res.lbLabelFontHeightF = 0.015              # default is a bit large

# Map options
res.mpLimitMode        = "Corners"
res.mpLeftCornerLatF   = lat[nlat-1][0] 
res.mpLeftCornerLonF   = lon[nlat-1][0]
res.mpRightCornerLatF  = lat[0][nlon-1]      
res.mpRightCornerLonF  = lon[0][nlon-1] 

res.mpGridAndLimbOn    = False
res.mpCenterLonF       = lon.Longitude_of_southern_pole 
res.mpCenterLatF       = lon.Latitude_of_southern_pole + 90
res.trYReverse         = True
res.tfDoNDCOverlay     = True

res.mpDataBaseVersion     = "MediumRes"
res.mpDataSetName         = "Earth..2"
res.mpOutlineBoundarySets = "AllBoundaries"

# Main Title
res.tiMainString      = filename

plot = Ngl.contour_map(wks,sbt,res)

Ngl.end()