; (11/19/2004 3.)(12/9/2004)(3/1/2005) ; ERA40 2D surface variables (ds120.0) conversion from GRIB to NetCDF ; ; NOTES ; ----- ; 1. These fields are on a regular 2.5 deg grid (144x73) ; 2. Field list: ; Var_name PN Name Long_name ; -------- -- ---- --------- ; SSHF_GDS0_SFC_113 146 SSHF Surface sensible heat flux ; SLHF_GDS0_SFC_113 147 SLHF Surface latent heat flux ; MSL_GDS0_SFC_123 151 MSL Mean sea level pressure ; TCC_GDS0_SFC_113 164 TCC Total cloud cover ; SSRD_GDS0_SFC_113 169 SSRD Surface solar radiation downwards ; STRD_GDS0_SFC_113 175 STRD Surface thermal radiation downwards ; SSR_GDS0_SFC_113 176 SSR Surface solar radiation ; STR_GDS0_SFC_113 177 STR Surface thermal radiation ; TSR_GDS0_SFC_113 178 TSR Top solar radiation ; TTR_GDS0_SFC_113 179 TTR Top thermal radiation ; EWSS_GDS0_SFC_113 180 EWSS East-West surface stress ; NSSS_GDS0_SFC_113 181 NSSS North-South surface stress ; LCC_GDS0_SFC_113 186 LCC Low cloud cover ; MCC_GDS0_SFC_113 187 MCC Medium cloud cover ; HCC_GDS0_SFC_113 188 HCC High cloud cover ; TSRC_GDS0_SFC_113 208 TSRC Top net solar radiation, clear sky ; TTRC_GDS0_SFC_113 209 TTRC Top net thermal radiation, clear sky ; SSRC_GDS0_SFC_113 210 SSRC Surface net solar radiation, clear sky ; STRC_GDS0_SFC_113 211 STRC Surface net thermal radiation, clear sky load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl" ; --------------- Subroutines --------------- undef ("grib_meta") procedure grib_meta(x) ; Bookkeeping: remove/add some GRIB meta-data begin if (isatt(x,"center")) then delete (x@center) end if if (isatt(x,"forecast_time")) then delete (x@forecast_time) end if if (isatt(x,"parameter_number")) then delete (x@parameter_number) end if if (isatt(x,"grid_number")) then delete (x@grid_number) end if if (isatt(x,"level_indicator")) then delete (x@level_indicator) end if if (isatt(x,"initial_time0")) then delete (x@initial_time0) end if if (isatt(x,"initial_time3")) then delete (x@initial_time3) end if x@_FillValue = 1.e20 x@missing_value = 1.e20 end ; --------------- Main block --------------- begin ; ..Read starting year, number of years and working directory from ; environment variable: setenv DS120 "${istrt};${nyrs};${wdir}" cla = getenv("DS120") clc = stringtochar(cla) lst = dimsizes(clc)-2 if (clc(lst) .eq. "/") then lst = lst-1 end if i = 0 do while (clc(i) .ne. ";") i = i+1 end do ie = i-1 i = i+1 do while (clc(i) .ne. ";") i = i+1 end do jb = ie+2 je = i-1 ke = je+2 istrt = stringtointeger(chartostring(clc(0:ie))) nfils = stringtointeger(chartostring(clc(jb:je))) iend = istrt+nfils-1 stp = False if (istrt .lt. 1958 .or. istrt .gt. 2001) then print("START YEAR: "+istrt+" Must have 1958 <= start year <= 2001") stp = True end if if (iend .lt. 1958 .or. iend .gt. 2001) then print("END YEAR : "+iend +" Must have 1958 <= end year <= 2001") stp = True end if if (stp) then exit end if tdir = chartostring(clc(ke:lst)) system ("mkdir -p "+tdir) wdir = tdir+"/" ; ..Create lists of MSS and local file names mssn = new((/nfils/),string) locn = new((/nfils/),string) idate = new((/nfils/),integer) i = -1 do iyr = istrt,istrt+nfils-1 i = i+1 mssn(i) = "/DSS/U0"+(iyr-629) locn(i) = "e4moda.sfc."+iyr idate(i) = iyr*10000+101 end do ; ..Read all files to local disk before proceeding with processing i = -1 do iyr = istrt,istrt+nfils-1 i = i+1 system ("msrcp -n mss:"+mssn(i)+" "+wdir+locn(i)) end do ; ..Regular grid nlat = 73 mlon = 144 lat = latGlobeF (nlat, "lat", "latitdue", "degrees_north") lon = lonGlobeF (mlon, "lon", "longitude", "degrees_east") ; ..Index of starting and ending times (0 --> 11 processes whole year) ntStrt = 0 ntLast = 11 ; ..Loop over all files do nf=0,nfils-1 print ("==========> nf="+nf+" <==========") wcStrt = systemfunc("date") ; wall clock start ; ..Open grib file grib_in = addfile(wdir+locn(nf)+".grb","r") wallClockElapseTime(wcStrt, "Read/open file: "+locn(nf), 0) ; ..Read and create times stime = grib_in->initial_time0(ntStrt:ntLast) time = grib_in->initial_time0_hours(ntStrt:ntLast) ntim = dimsizes(time) DATE = doubletointeger(grib_in->initial_time0_encoded(ntStrt:ntLast)) date = DATE/100 date@long_name = "current date as 8 digit integer (YYYYMMDD)" date!0= "time" print(date) if (date(ntStrt) .eq. idate(nf)) then print("Date in input file matches specified date") else print("STOP -- date in input file does not match specified date") exit() end if hh = (DATE-date*100) datesec = hh*3600 datesec!0 = "time" datesec@long_name = "current seconds of current date" datesec@units = "s" if (nf.eq.0) then ; 1st file only [debug] print (stime) print (time) print (DATE) print (date) print (datesec) end if ; ..Read standard scalar variables: GRIB lats run N->S SSHF = grib_in->SSHF_GDS0_SFC_113(ntStrt:ntLast,:,:) SLHF = grib_in->SLHF_GDS0_SFC_113(ntStrt:ntLast,:,:) MSL = grib_in->MSL_GDS0_SFC_123(ntStrt:ntLast,:,:) TCC = grib_in->TCC_GDS0_SFC_113(ntStrt:ntLast,:,:) SSRD = grib_in->SSRD_GDS0_SFC_113(ntStrt:ntLast,:,:) STRD = grib_in->STRD_GDS0_SFC_113(ntStrt:ntLast,:,:) SSR = grib_in->SSR_GDS0_SFC_113(ntStrt:ntLast,:,:) STR = grib_in->STR_GDS0_SFC_113(ntStrt:ntLast,:,:) TSR = grib_in->TSR_GDS0_SFC_113(ntStrt:ntLast,:,:) TTR = grib_in->TTR_GDS0_SFC_113(ntStrt:ntLast,:,:) EWSS = grib_in->EWSS_GDS0_SFC_113(ntStrt:ntLast,:,:) NSSS = grib_in->NSSS_GDS0_SFC_113(ntStrt:ntLast,:,:) LCC = grib_in->LCC_GDS0_SFC_113(ntStrt:ntLast,:,:) MCC = grib_in->MCC_GDS0_SFC_113(ntStrt:ntLast,:,:) HCC = grib_in->HCC_GDS0_SFC_113(ntStrt:ntLast,:,:) TSRC = grib_in->TSRC_GDS0_SFC_113(ntStrt:ntLast,:,:) TTRC = grib_in->TTRC_GDS0_SFC_113(ntStrt:ntLast,:,:) SSRC = grib_in->SSRC_GDS0_SFC_113(ntStrt:ntLast,:,:) STRC = grib_in->STRC_GDS0_SFC_113(ntStrt:ntLast,:,:) ; ..Make lat coords S->N SSHF = SSHF(:,::-1,:) SLHF = SLHF(:,::-1,:) MSL = MSL(:,::-1,:) TCC = TCC(:,::-1,:) SSRD = SSRD(:,::-1,:) STRD = STRD(:,::-1,:) SSR = SSR(:,::-1,:) STR = STR(:,::-1,:) TSR = TSR(:,::-1,:) TTR = TTR(:,::-1,:) EWSS = EWSS(:,::-1,:) NSSS = NSSS(:,::-1,:) LCC = LCC(:,::-1,:) MCC = MCC(:,::-1,:) HCC = HCC(:,::-1,:) TSRC = TSRC(:,::-1,:) TTRC = TTRC(:,::-1,:) SSRC = SSRC(:,::-1,:) STRC = STRC(:,::-1,:) ; ..Add CAM2-like meta-data SSHF!0 = "time" SLHF!0 = "time" MSL!0 = "time" TCC!0 = "time" SSRD!0 = "time" STRD!0 = "time" SSR!0 = "time" STR!0 = "time" TSR!0 = "time" TTR!0 = "time" EWSS!0 = "time" NSSS!0 = "time" LCC!0 = "time" MCC!0 = "time" HCC!0 = "time" TSRC!0 = "time" TTRC!0 = "time" SSRC!0 = "time" STRC!0 = "time" SSHF!1 = "lat" SLHF!1 = "lat" MSL!1 = "lat" TCC!1 = "lat" SSRD!1 = "lat" STRD!1 = "lat" SSR!1 = "lat" STR!1 = "lat" TSR!1 = "lat" TTR!1 = "lat" EWSS!1 = "lat" NSSS!1 = "lat" LCC!1 = "lat" MCC!1 = "lat" HCC!1 = "lat" TSRC!1 = "lat" TTRC!1 = "lat" SSRC!1 = "lat" STRC!1 = "lat" SSHF!2 = "lon" SLHF!2 = "lon" MSL!2 = "lon" TCC!2 = "lon" SSRD!2 = "lon" STRD!2 = "lon" SSR!2 = "lon" STR!2 = "lon" TSR!2 = "lon" TTR!2 = "lon" EWSS!2 = "lon" NSSS!2 = "lon" LCC!2 = "lon" MCC!2 = "lon" HCC!2 = "lon" TSRC!2 = "lon" TTRC!2 = "lon" SSRC!2 = "lon" STRC!2 = "lon" ; ..Remove GRIB meta-data grib_meta(SSHF) grib_meta(SLHF) grib_meta(MSL) grib_meta(TCC) grib_meta(SSRD) grib_meta(STRD) grib_meta(SSR) grib_meta(STR) grib_meta(TSR) grib_meta(TTR) grib_meta(EWSS) grib_meta(NSSS) grib_meta(LCC) grib_meta(MCC) grib_meta(HCC) grib_meta(TSRC) grib_meta(TTRC) grib_meta(SSRC) grib_meta(STRC) ; ..Set fill values ; SSHF@_FillValue = 1.e20 ; SLHF@_FillValue = 1.e20 ; MSL@_FillValue = 1.e20 ; TCC@_FillValue = 1.e20 ; SSRD@_FillValue = 1.e20 ; STRD@_FillValue = 1.e20 ; SSR@_FillValue = 1.e20 ; STR@_FillValue = 1.e20 ; TSR@_FillValue = 1.e20 ; TTR@_FillValue = 1.e20 ; EWSS@_FillValue = 1.e20 ; NSSS@_FillValue = 1.e20 ; LCC@_FillValue = 1.e20 ; MCC@_FillValue = 1.e20 ; HCC@_FillValue = 1.e20 ; TSRC@_FillValue = 1.e20 ; TTRC@_FillValue = 1.e20 ; SSRC@_FillValue = 1.e20 ; STRC@_FillValue = 1.e20 ; ..Print out summary of variables printVarSummary(SSHF) printVarSummary(SLHF) printVarSummary(MSL) printVarSummary(TCC) printVarSummary(SSRD) printVarSummary(STRD) printVarSummary(SSR) printVarSummary(STR) printVarSummary(TSR) printVarSummary(TTR) printVarSummary(EWSS) printVarSummary(NSSS) printVarSummary(LCC) printVarSummary(MCC) printVarSummary(HCC) printVarSummary(TSRC) printVarSummary(TTRC) printVarSummary(SSRC) printVarSummary(STRC) ; ..Open output NetCDF file; create parameters and output file system ("rm -f "+wdir+locn(nf)+".nc") netcdf_out = addfile(wdir+locn(nf)+".nc","c") ; ..Assign file attributes wcDef = systemfunc("date") nl = inttochar(10) ; carriage return fAtt = True fAtt@title = "ERA40 GRIB-to-NetCDF" fAtt@story = nl+\ "ERA40 GRIB file: These are monthly mean 2D surface variables on a 2.5 "+nl+\ "degree regular grid (144x73). We simply read in the 19 variables that "+nl+\ "we want and write them out to another NetCDF file after changing or "+nl+\ "removing some of the meta-data. "+nl fAtt@Conventions = "None" fAtt@source_file = locn(nf) fAtt@source_center = "European Center for Medium-Range Weather Forecasts - Reading" fAtt@creation_date = systemfunc("date") fileattdef( netcdf_out, fAtt ) ; ..Predefine coordinate information dimNames = (/ "time", "lat", "lon" /) dimSizes = (/ ntim , nlat , mlon /) dimUnlim = (/ True , False, False /) filedimdef(netcdf_out, dimNames, dimSizes, dimUnlim) filevardef (netcdf_out, "time", typeof(time), "time") filevarattdef(netcdf_out, "time", time) filevardef (netcdf_out, "lat" , typeof(lat) , "lat") filevarattdef(netcdf_out, "lat" , lat) filevardef (netcdf_out, "lon" , typeof(lon) , "lon") filevarattdef(netcdf_out, "lon" , lon) ; ..Predefine variable sizes filevardef (netcdf_out, "date" , typeof(date) , getvardims(date)) filevarattdef(netcdf_out, "date" , date) filevardef (netcdf_out, "datesec", typeof(datesec), getvardims(datesec)) filevarattdef(netcdf_out, "datesec", datesec) filevardef (netcdf_out, "SSHF" , typeof(SSHF) , getvardims(SSHF)) filevarattdef(netcdf_out, "SSHF" , SSHF) filevardef (netcdf_out, "SLHF" , typeof(SLHF) , getvardims(SLHF)) filevarattdef(netcdf_out, "SLHF" , SLHF) filevardef (netcdf_out, "MSL" , typeof(MSL) , getvardims(MSL)) filevarattdef(netcdf_out, "MSL" , MSL) filevardef (netcdf_out, "TCC" , typeof(TCC) , getvardims(TCC)) filevarattdef(netcdf_out, "TCC" , TCC) filevardef (netcdf_out, "SSRD" , typeof(SSRD) , getvardims(SSRD)) filevarattdef(netcdf_out, "SSRD" , SSRD) filevardef (netcdf_out, "STRD" , typeof(STRD) , getvardims(STRD)) filevarattdef(netcdf_out, "STRD" , STRD) filevardef (netcdf_out, "SSR" , typeof(SSR) , getvardims(SSR)) filevarattdef(netcdf_out, "SSR" , SSR) filevardef (netcdf_out, "STR" , typeof(STR) , getvardims(STR)) filevarattdef(netcdf_out, "STR" , STR) filevardef (netcdf_out, "TSR" , typeof(TSR) , getvardims(TSR)) filevarattdef(netcdf_out, "TSR" , TSR) filevardef (netcdf_out, "TTR" , typeof(TTR) , getvardims(TTR)) filevarattdef(netcdf_out, "TTR" , TTR) filevardef (netcdf_out, "EWSS" , typeof(EWSS) , getvardims(EWSS)) filevarattdef(netcdf_out, "EWSS" , EWSS) filevardef (netcdf_out, "NSSS" , typeof(NSSS) , getvardims(NSSS)) filevarattdef(netcdf_out, "NSSS" , NSSS) filevardef (netcdf_out, "LCC" , typeof(LCC) , getvardims(LCC)) filevarattdef(netcdf_out, "LCC" , LCC) filevardef (netcdf_out, "MCC" , typeof(MCC) , getvardims(MCC)) filevarattdef(netcdf_out, "MCC" , MCC) filevardef (netcdf_out, "HCC" , typeof(HCC) , getvardims(HCC)) filevarattdef(netcdf_out, "HCC" , HCC) filevardef (netcdf_out, "TSRC" , typeof(TSRC) , getvardims(TSRC)) filevarattdef(netcdf_out, "TSRC" , TSRC) filevardef (netcdf_out, "TTRC" , typeof(TTRC) , getvardims(TTRC)) filevarattdef(netcdf_out, "TTRC" , TTRC) filevardef (netcdf_out, "SSRC" , typeof(SSRC) , getvardims(SSRC)) filevarattdef(netcdf_out, "SSRC" , SSRC) filevardef (netcdf_out, "STRC" , typeof(STRC) , getvardims(STRC)) filevarattdef(netcdf_out, "STRC" , STRC) wallClockElapseTime(wcDef, "Definition phase took ", 0) ; ..Write data values to predefined locations wcWrite = systemfunc("date") netcdf_out->time = (/ time /) netcdf_out->lat = (/ lat /) netcdf_out->lon = (/ lon /) netcdf_out->date = (/ date /) netcdf_out->datesec = (/ datesec /) netcdf_out->SSHF = (/ SSHF /) netcdf_out->SLHF = (/ SLHF /) netcdf_out->MSL = (/ MSL /) netcdf_out->TCC = (/ TCC /) netcdf_out->SSRD = (/ SSRD /) netcdf_out->STRD = (/ STRD /) netcdf_out->SSR = (/ SSR /) netcdf_out->STR = (/ STR /) netcdf_out->TSR = (/ TSR /) netcdf_out->TTR = (/ TTR /) netcdf_out->EWSS = (/ EWSS /) netcdf_out->NSSS = (/ NSSS /) netcdf_out->LCC = (/ LCC /) netcdf_out->MCC = (/ MCC /) netcdf_out->HCC = (/ HCC /) netcdf_out->TSRC = (/ TSRC /) netcdf_out->TTRC = (/ TTRC /) netcdf_out->SSRC = (/ SSRC /) netcdf_out->STRC = (/ STRC /) wallClockElapseTime(wcWrite, "Write phase took ", 0) delete (stime) ; may change in next file delete (time) delete (date) delete (datesec) delete (SSHF) delete (SLHF) delete (MSL) delete (TCC) delete (SSRD) delete (STRD) delete (SSR) delete (STR) delete (TSR) delete (TTR) delete (EWSS) delete (NSSS) delete (LCC) delete (MCC) delete (HCC) delete (TSRC) delete (TTRC) delete (SSRC) delete (STRC) end do ; loop on files print("Normal termination of program ds120.0.ncl") end