load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl" load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl" load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl" load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl" begin ;## path to wpp outfiles path = "/import/wrkdir1/asemenov/ERA_INTERIM_NEW/NOAA_optimum_interpolation_analysis/15_March/" ;## file name files = systemfunc("ls "+path+"avhrr-only-v2.*") ;## create array s = new ( (/30,1,1,720,1440/), "float" ) s@_FillValue = -99999 s = s@_FillValue s@long_name = "sst" s@units = "oC" s!0 = "file_n" s!1 = "time" s!2 = "zlev" s!3 = "lat" s!4 = "lon" printVarSummary(s) i = new ( (/30,1,1,720,1440/), "float" ) i@_FillValue = -99999 i = i@_FillValue i@long_name = "sea ice concentration" i@units = "%" i!0 = "file_n" i!1 = "time" i!2 = "zlev" i!3 = "lat" i!4 = "lon" printVarSummary(i) ;## read throught the files do j = 0, dimsizes(files)-1 fname = files(j) ;## open file print("Opening file: "+fname+".nc") f = addfile (fname+".nc", "r") lat=f->lat lon=f->lon zlev=f->zlev time=f->time ; sa=f->sst sa=short2flt(f->sst) printVarSummary(sa) ; ic=f->ice ic=short2flt(f->ice) printVarSummary(ic) s(j,:,:,:,:)=sa(:,:,:,:) i(j,:,:,:,:)=ic(:,:,:,:) end do ;printVarSummary(s) ;printMinMax(s,True) ;printVarSummary(i) ;printMinMax(i,True) S_avg=s I_avg=i do t =0,719 if (lat(t).gt.56 .and. lat(t).le.57.5) then I_avg(29,:,:,t,:)=0.1 S_avg(29,:,:,t,:)=0.9 end if if (lat(t).gt.57.5 .and. lat(t).le.59.) then I_avg(29,:,:,t,:)=0.2 S_avg(29,:,:,t,:)=0.7 end if if (lat(t).gt.59. .and. lat(t).le.59.5) then I_avg(29,:,:,t,:)=0.3 S_avg(29,:,:,t,:)=0.6 end if if (lat(t).gt.59.5 .and. lat(t).le.60.150) then I_avg(29,:,:,t,:)=0.5 S_avg(29,:,:,t,:)=0.5 end if if (lat(t).gt.60.150 .and. lat(t).le.60.85) then I_avg(29,:,:,t,:)=0.7 S_avg(29,:,:,t,:)=0.3 end if if (lat(t).gt.60.85 .and. lat(t).le.61.4) then I_avg(29,:,:,t,:)=0.8 S_avg(29,:,:,t,:)=0.2 end if if (lat(t).gt.61.4 .and. lat(t).le.61.9) then I_avg(29,:,:,t,:)=0.9 S_avg(29,:,:,t,:)=0.1 end if if (lat(t).gt.61.9) then I_avg(29,:,:,t,:)=1.0 S_avg(29,:,:,t,:)=0.0 end if end do ;S_avg=dim_avg_n_Wrap(s,0) ;I_avg=dim_avg_n_Wrap(i,0) ;printVarSummary(s) ;printVarsummary(S_avg) ;rintVarSummary(S_avg) ;printMinMax(I_avg,True) printVarSummary(I_avg) printMinMax(I_avg,True) wks = gsn_open_wks("pdf" ,"polar_1982") ; open a ps file ;gsn_define_colormap(wks,"BlueWhiteOrangeRed") ;gsn_merge_colormaps(wks,"precip3_16lev","WhiteYellowOrangeRed") gsn_define_colormap(wks,"WhViBlGrYeOrRe") res = True ; plot mods desired res@gsnPolar = "NH" ; specify the hemisphere res@mpMinLatF = 55 ; minimum lat to plot res@mpFillOn = False res@cnFillOn = True ; color fill res@cnLinesOn = False ;res@cnLevelSelectionMode = "ManualLevels" ; res@cnMinLevelValF = 0. ; res@cnMaxLevelValF = 1. ; res@cnLevelSpacingF = 0.1 ; res@gsnSpreadColors = False ; delete (cres@cnFillColors) ; res@cnFillColors = (/ 6,4,75,91,107,123,139,155,171,187,193,208,218,225,234 /) res@cnLevelSelectionMode = "ManualLevels" res@cnMinLevelValF = 0. res@cnMaxLevelValF = 1. res@cnLevelSpacingF = 0.1 res@gsnSpreadColors = False ;res@gsnMaximize = True ; delete (cres@cnFillColors) res@cnFillColors = (/ 2,4,6,8,10,12,14,15,16,18,20,22 /) res@mpGridAndLimbOn = True res@mpGridLineDashPattern = 2 res@mpGeophysicalLineColor = "Black" res@mpGeophysicalLineThicknessF = 1.5 res@mpNationalLineColor = "Black" res@mpNationalLineThicknessF = 1.5 res@mpUSStateLineColor = "Black" res@mpUSStateLineThicknessF = 1.5 res@mpPerimLineColor = "Black" res@mpDataBaseVersion = "MediumRes" res@mpFillOn = True res@mpFillColors = (/-1,-1,0,-1/) res@mpGridLatSpacingF = 30. res@mpGridLonSpacingF = 30. res@mpOutlineOn = True ; res@pmTickMarkDisplayMode = "Always" ; res@tmXTOn = False ; res@tmXBOn = True ; res@tmYLOn = False ; res@tmYROn = False ; res@tmXTLabelFontHeightF = 0.025 ; res@tmXBLabelFontHeightF = 0.025 ; res@tmYLLabelFontHeightF = 0.025 ; res@tmYRLabelFontHeightF = 0.025 ; res@tmYLValues = (/ 60,70,80 /) ; res@tmXBValues = (/ 150., -180., -150. /) ; res@tmXBLabelDeltaF = -0.8 ; res@tmYLLabelDeltaF = -0.8 ;LABEL BAR cresOURCES ; res@pmLabelBarDisplayMode = "Always" ; res@lbLabelAutoStride = True ; res@lbOrientation = "Vertical" res@lbLabelFontHeightF = 0.025 res@lbPerimOn = False res@lbLabelFontHeightF = 0.025 ; res@lbAutoManage = False ; res@cnInfoLabelOn = False ;GSN cresOURCES res@gsnLeftString = " " res@gsnRightString = " " res@gsnMaximize = True ; res@gsnDraw = False ; res@gsnFrame = False ; ;TITLE cresOURCES res@tiMainOn = True res@tiMainPosition = "Center" res@tiMainFontHeightF = 0.025 a = addfile("$NCARG_ROOT/lib/ncarg/data/cdf/landsea.nc","r") lsm_I_avg = landsea_mask(a->LSMASK,I_avg&lat,I_avg&lon) ; read in land sea mask, and pass it and ; the t85 lat/lon arrays into landsea_mask I_avg = mask(I_avg,lsm_I_avg.eq.1,False) lsm_S_avg = landsea_mask(a->LSMASK,S_avg&lat,S_avg&lon) ; read in land sea mask, and pass it and ; the t85 lat/lon arrays into landsea_mask S_avg = mask(S_avg,lsm_S_avg.eq.1,False) ;ata = mask(data,lsm.ge.1,False) res@tiMainString = "SST for extended sea ice experiment,~C~March 15 2011, C" ; plot = gsn_csm_contour_map_polar(wks,s(0,0,0,:,:),res) ; ; plot = gsn_csm_contour_map_polar(wks,S_avg(29,0,0,:,:),res) S_avg=S_avg+273 ;=================================================================== ; Assume variables T, PS and ORO exist and that they have ; associated meta data: (a) coordinate variables time, lev, lat, lon ; and (b) attributes ;=================================================================== ntim = dimsizes(time) ; get dimension sizes nlev = dimsizes(zlev) nlat = dimsizes(lat) nlon = dimsizes(lon) diro = "./" ; Output directory filo = "double.nc" ; Output file system("/bin/rm -f " + diro + filo) ; remove if exists fout = addfile (diro + filo, "c") ; open output file ;=================================================================== ; explicitly declare file definition mode. Improve efficiency. ;=================================================================== setfileoption(fout,"DefineMode",True) ;=================================================================== ; create global attributes of the file ;=================================================================== fAtt = True ; assign file attributes fAtt@title = "NCL Efficient Approach to netCDF Creation" fAtt@source_file = "avhrr-only-v2.20110315.nc" fAtt@Conventions = "None" fAtt@creation_date = systemfunc ("date") fileattdef( fout, fAtt ) ; copy file attributes ;=================================================================== ; predefine the coordinate variables and their dimensionality ; Note: to get an UNLIMITED record dimension, we set the dimensionality ; to -1 (or the actual size) and set the dimension name to True. ;=================================================================== dimNames = (/"time", "zlev", "lat", "lon"/) dimSizes = (/ -1 , nlev, nlat, nlon /) dimUnlim = (/ True , False, False, False/) filedimdef(fout,dimNames,dimSizes,dimUnlim) ;=================================================================== ; predefine the the dimensionality of the variables to be written out ;=================================================================== ; Here we are using NCL functions to facilitate defining ; each variable's dimension name(s) and type. ; The following could be replaced with explicit, user defined dimension ; names different from those associated with the variable in memory. ; Say, PS(time,lat,lon) in the NCL script. They could be redefined for the file via: ; filevardef(fout, "PS" ,typeof(PS) ,(/"TIME","latitude","longitude"/)) ;=================================================================== filevardef(fout, "time" ,typeof(time),getvardims(time)) filevardef(fout, "zlev" ,typeof(zlev),getvardims(zlev) ) filevardef(fout, "lat" ,typeof(lat),getvardims(lat)) filevardef(fout, "lon" ,typeof(lon),getvardims(lon)) filevardef(fout, "S_avg" ,typeof(S_avg) , (/"time","zlev","lat","lon"/) ) filevardef(fout, "I_avg" ,typeof(I_avg) , (/"time","zlev","lat","lon"/) ) ; filevardef(fout, "TOPOG",typeof(ORO),getvardims(ORO)) ; variable name on the file ; different from name on script ;=================================================================== ; Copy attributes associated with each variable to the file ; All attributes associated with each variable will be copied. ;==================================================================== filevarattdef(fout,"S_avg",S_avg) ; copy T attributes filevarattdef(fout,"time" ,time) ; copy time attributes filevarattdef(fout,"zlev" ,zlev) ; copy lev attributes filevarattdef(fout,"lat" ,lat) ; copy lat attributes filevarattdef(fout,"lon" ,lon) ; copy lon attributes filevarattdef(fout,"I_avg" ,I_avg) ; copy PS attributes ; filevarattdef(fout,"TOPOG",ORO) ; copy TOPOG attributes ;=================================================================== ; explicitly exit file definition mode. **NOT REQUIRED** ;=================================================================== setfileoption(fout,"DefineMode",False) ;=================================================================== ; output only the data values since the dimensionality and such have ; been predefined. The "(/", "/)" syntax tells NCL to only output the ; data values to the predefined locations on the file. ;==================================================================== fout->time = (/time/) fout->zlev = (/zlev/) fout->lat = (/lat/) fout->lon = (/lon/) fout->I_avg = (/I_avg/) fout->S_avg = (/S_avg/) ; fout->TOPOG = (/ORO/) end