# # File: # hdf1.py # # Synopsis: # Plots EOS-DIS data. # # Category: # Contours over maps # Labelbar # Maps # # Author: # Mary Haley (based on NCL example from Dennis Shea) # # Date of initial publication: # July, 2011 # # Description: # This example uses Ngl.dim_gbits to extra bit information from # cloud mask data. # # Effects illustrated: # o Reading from a NetCDF file. # o Setting map resources. # o Setting labelbar resources. # o Contouring in triangular mesh mode. # # Output: # A single visualization is produced showing cloud mask data # over a satellite projection. # # #****************************************************************** # Bit fields within each byte are numbered from the left: # 7, 6, 5, 4, 3, 2, 1, 0. # The left-most bit (bit 7) is the most significant bit. # The right-most bit (bit 0) is the least significant bit. # # bit field Description Key # --------- ----------- --- # 2, 1 Unobstructed FOV Quality Flag 00 (0) = Cloudy # 01 (1) = Uncertain # 10 (2) = Probably Clear # 11 (3) = Clear #---Import support modules from __future__ import print_function import numpy, os, sys, Nio, Ngl #---Test if file exists flnm = "MOD05_L2.A2002161.1830.004.2003221153410.hdf" if(not os.path.exists(flnm)): print("You do not have the necessary {} file to run this example.".format(flnm)) print("You can get the file from") print(" http://www.ncl.ucar.edu/Applications/Data/#hdf") sys.exit() #---Open netCDF file hfile = Nio.open_file(flnm) #---Read lat/lon [5km] and Cloud_Mask [1km] cldmsk = hfile.variables["Cloud_Mask_QA"][:] lat2d = hfile.variables["Latitude"][:] lon2d = hfile.variables["Longitude"][:] dimll = lat2d.shape nlat = dimll[0] # "along" swath [5 km] nlon = dimll[1] # "across" swath [5 km] # # Latitude/Longitude arrays are every 5km. Cloud_Mask is every 1km # It is necessary to decimate the Cloud_Mask to match lat/lon # dimcm = cldmsk.shape NY = dimcm[0] # "along" swath [1 km] MX = dimcm[1] # "across" swath [1 km] cm = cldmsk[4::5,4::5] # decimate to match lat2d/lon2d dimsizes del cldmsk # no longer needed # # Extract the desired bits and attach lat/lon arrays # cmfov = 00 [0] cloudy 01 [1] uncertain # 10 [2] probably clear 00 [3] confidently clear # tmp = Ngl.dim_gbits(numpy.ravel(cm),5,2,6,nlat*nlon) cmfov = numpy.reshape(tmp, [nlat,nlon]) del tmp del cm #---Open a PNG file wks_type = "png" wks = Ngl.open_wks(wks_type,"hdf1") #---Create variable to hold list of plot resources. res = Ngl.Resources() colors = ["gray","red","green","white"] res.nglDraw = False # Don't create plot or res.nglFrame = False # advance frame yet. res.cnFillOn = True # turn on color fill res.cnLinesOn = False # turn off contour lines res.cnLineLabelsOn = False # turn of contour labels res.cnFillMode = "RasterFill" # turn on raster mode res.cnFillPalette = colors res.pmLabelBarDisplayMode = "Never" # Turn off labelbar res.cnLevelSelectionMode = "ManualLevels" # set manual contour levels res.cnMinLevelValF = 1 # set min contour level res.cnMaxLevelValF = 3 # one less than max res.cnLevelSpacingF = 1 # set contour spacing res.mpDataBaseVersion = "MediumRes" res.mpOutlineBoundarySets= "National" res.mpProjection = "Satellite" # choose map projection res.mpFillOn = False # turn off map fill res.mpGridAndLimbOn = False # turn off lat/lon lines res.mpCenterLonF = (numpy.min(lon2d)+numpy.max(lon2d))*0.5 res.mpCenterLatF = (numpy.min(lat2d)+numpy.max(lat2d))*0.5 res.mpLimitMode = "LatLon" res.mpMinLatF = numpy.min(lat2d)-1 # min lat res.mpMaxLatF = numpy.max(lat2d)+1 # max lat res.mpMinLonF = numpy.min(lon2d)-1 # min lon res.mpMaxLonF = numpy.max(lon2d)+1 # max lon res.tiMainString = "Unobstructed FOV Quality Flag" res.trGridType = "TriangularMesh" # faster graphic rendering res.sfXArray = lon2d # 5 km 2D lat/lon arrays res.sfYArray = lat2d #---Create plot (it won't get drawn yet) plot = Ngl.contour_map(wks,cmfov,res) #---Resources for creating a labelbar lbres = Ngl.Resources() lbres.nglDraw = False lbres.lbPerimOn = False # turn off perimeter. lbres.vpWidthF = 0.60 # width and height lbres.vpHeightF = 0.15 lbres.lbOrientation = "Horizontal" lbres.lbLabelPosition = "Center" # label position lbres.lbFillColors = colors lbres.lbMonoFillPattern = True # solid color fill lbres.lbLabelFontHeightF = 0.02 lbres.lbTitleOn = True lbres.lbTitleString = "0=cldy, 1=uncertain, 2=prob clr, 3=clr" lbres.lbTitlePosition = "Bottom" lbres.lbTitleFontHeightF = 0.02 lbres.lbTitleOffsetF = 0.00 #---Create the labelbar labels = ["0","1","2","3"] lbid = Ngl.labelbar_ndc (wks,len(labels),labels,0.0,0.,lbres) # # Create some annotation resources indicating how we want to # attach the labelbar to the plot. The default is the center # of the plot. Below amOrthogonalPosF is set to move the # labelbar down and outside the plot. # amres = Ngl.Resources() amres.amOrthogonalPosF = 0.7 annoid = Ngl.add_annotation(plot,lbid,amres) txres = Ngl.Resources() txres.txFontHeightF = 0.01 txres.txJust = "TopRight" txres.txPerimOn = True txid = Ngl.add_text(wks,plot,flnm,res.mpMaxLonF,res.mpMaxLatF,txres) #---This will resize plot so it and the labelbar fit in the frame. Ngl.maximize_plot(wks, plot) #---Drawing the original map also draws the attached labelbar. Ngl.draw(plot) Ngl.frame(wks) Ngl.end()