; Program to examine the HAMSR thinned data on 5x5 km grid ; Andrew Kren ; May 21, 2018 begin path_hamsr = "/scratch4/BMC/shout/Andrew.Kren/HAMSR/" ; location of HAMSR data ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; end of user's settings ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; fils = systemfunc("ls " + path_hamsr + "*.nc") ; file paths w = addfiles(fils,"r") ; read data lat = w[:]->lat lat := lat * 0.001 lon = w[:]->lon lon := lon *0.001 time = w[:]->time t = short2flt(w[:]->ham_airT) q = short2flt(w[:]->ham_airQ) rh = short2flt(w[:]->ham_airRH) lev = short2flt(w[:]->ham_pres_levels) lev := lev(0:24) t@units = "K" ; calculate dates time_hamsr = cd_calendar(time,0) year = time_hamsr(:,0) month = time_hamsr(:,1) day = time_hamsr(:,2) hour = time_hamsr(:,3) minute = time_hamsr(:,4) sec = time_hamsr(:,5) ; calculate specific humidity shum = mixhum_ptrh(conform(t,lev,1),t,rh,-2) ; g/kg ; bin the data by level do i=0,dimsizes(lev)-1 ; loop over levls wks = gsn_open_wks("png","hamsr_t_histo_" + lev(i) + "mb") ; send graphics to PNG file gsn_define_colormap(wks,"temp1") ; choose colormap res = True res@tiMainString = "HAMSR Temperature at " + lev(i) + " hPa" min_val = toint(floor(min(t(:,i)))) max_val = toint(ceil(max(t(:,i)))) res@gsnHistogramClassIntervals = ispan(min_val,max_val,2) res@gsnHistogramComputePercentages = True ; change left axis to % res@tiXAxisString = "Temperature (K)" res@tmXBLabelAngleF = 315. ; change label angle plot=gsn_histogram(wks,t(:,i),res) ; create histogram with 10 bins delete([/res,min_val,max_val/]) wks = gsn_open_wks("png","hamsr_q_histo_" + lev(i) + "mb") ; send graphics to PNG file gsn_define_colormap(wks,"temp1") ; choose colormap res = True res@tiMainString = "HAMSR Specific Humidity at " + lev(i) + " hPa" min_val = min(shum(:,i)) max_val = max(shum(:,i)) temp = fspan(min_val,max_val,10) res@gsnHistogramClassIntervals = (/sprintf("%4.1f",temp)/) res@gsnHistogramComputePercentages = True ; change left axis to % res@tiXAxisString = "Specific Humidity (g/kg)" res@tmXBLabelAngleF = 315. ; change label angle plot=gsn_histogram(wks,shum(:,i),res) ; create histogram with 10 bins delete([/res,min_val,max_val,temp/]) wks = gsn_open_wks("png","hamsr_rh_histo_" + lev(i) + "mb") ; send graphics to PNG file gsn_define_colormap(wks,"temp1") ; choose colormap res = True res@tiMainString = "HAMSR Relative Humidity at " + lev(i) + " hPa" min_val = toint(floor(min(rh(:,i)))) max_val = toint(ceil(max(rh(:,i)))) if (mod(min_val,10).ge.5) then min_val = min_val / 10 * 10 + 10 else min_val = min_val / 10 * 10 end if if (mod(max_val,10).ge.5) then max_val = max_val / 10 * 10 + 10 else max_val = max_val / 10 * 10 end if res@gsnHistogramClassIntervals = ispan(min_val,max_val,5) res@gsnHistogramComputePercentages = True ; change left axis to % res@tmXBLabelAngleF = 315. ; change label angle res@tiXAxisString = "Relative Humidity (%)" plot=gsn_histogram(wks,rh(:,i),res) ; create histogram with 10 bins delete([/res,min_val,max_val/]) end do ; plot contour of data on x-y plot t_reorder = t(HAMSR_levels|:,time|:) t_reorder!0 = "lev" t_reorder&lev = lev shum!0 = "time" shum!1 = "lev" shum&lev = lev shum_reorder = shum(lev|:,time|:) shum_reorder@long_name = "HAMSR Specific Humidity" shum_reorder@units = "g/kg" wks = gsn_open_wks ("png", "hamsr_t_contour") ; send graphics to PNG file res = True ; plot mods desired res@cnFillOn = True ; turn on color res@cnFillPalette = "BlAqGrYeOrReVi200" ; set color map res@cnLinesOn = False ; no contour lines res@cnLevelSelectionMode = "ManualLevels" ; set manual contour levels res@cnMinLevelValF = 240.0 ; set min contour level res@cnMaxLevelValF = 300.0 ; set max contour level res@cnLevelSpacingF = 2.0 ; set contour spacing res@tiYAxisString = "Pressure (hPa)" res@tiXAxisString = "Time" plot = gsn_csm_pres_hgt (wks,t_reorder(0:20,0:100),res) wks = gsn_open_wks ("png", "hamsr_q_contour") ; send graphics to PNG file res = True ; plot mods desired res@cnFillOn = True ; turn on color res@cnFillPalette = "BlAqGrYeOrReVi200" ; set color map res@cnLinesOn = False ; no contour lines res@cnLevelSelectionMode = "ManualLevels" ; set manual contour levels res@cnMinLevelValF = 0.5 ; set min contour level res@cnMaxLevelValF = 22.0 ; set max contour level res@cnLevelSpacingF = 0.5 ; set contour spacing res@tiYAxisString = "Pressure (hPa)" res@tiXAxisString = "Time" plot = gsn_csm_pres_hgt (wks,shum_reorder(0:20,0:100),res) delete(res) ; plot data spatially outfile = "hamsr_t_spatial_map" + "_" + tostring(ymd(0)) wks = gsn_open_wks("png",outfile) gsn_define_colormap(wks,"BlueWhiteOrangeRed") res = True res@pmTickMarkDisplayMode = "Always" ; turn on fancy tickmarks res@mpOutlineOn = True ; turn the map outline on res@mpDataBaseVersion = "MediumRes" ; Medium resolution res@mpDataSetName = "Earth..4" res@gsnMaximize = True res@tiMainString = "HAMSR Temperature from " + tostring(year(0)) + tostring(month(0)) + tostring(day(0)) + tostring(hour(0)) + tostring(minute(0)) + "." + tostring(sec(0)) + " to " + tostring(year(100)) + tostring(month(100)) + tostring(day(100)) + tostring(hour(100)) + tostring(minute(100)) + "." + tostring(sec(100)) ;---Zoom in on United States. res@mpMinLatF = 20 res@mpMaxLatF = 45 res@mpMinLonF = 270 res@mpMaxLonF = 300 res@mpCenterLonF = 285 res@mpFillOn = False res@gsnRightString = "" res@gsnLeftString = "" res@mpGridAndLimbOn = True res@mpGridLatSpacingF = 10.0 res@mpGridLonSpacingF = 10.0 res@mpOutlineBoundarySets = "GeophysicalAndUSStates" res@mpGeophysicalLineThicknessF = 2.5 ; thickness of boundaries res@mpUSStateLineThicknessF = 2.5 res@gsnFrame = False ; don't advance frame res@gsnDraw = False ; don't draw plot plot = gsn_csm_map(wks,res) gres = True gres@gsLineColor = "Red" gres@gsLineThicknessF = 3.0 dum3 = gsn_add_polyline(wks,plot,lon(0:100),lat(0:100),gres) draw(plot) frame(wks) end