Dear NCL community,
First, I want to thank MARY for her suggestion. I have tried Mary's ways ,but I failed in outputing data. And I don't know how to deal with some errors(WARNING).
Scripts and results below are in MARY's way and in an other way with a loop . Kindly hope someone give me more suggestions.
Thanks for your attention.
Echo
NO.1(Mary's way):
begin
;---Input file
diri = "./"
fili = "ACCESS1-0/Historial/r1i1p1/tas_Amon_ACCESS1-0_historical_r1i1p1_185001-200512.nc"
year = 1901 ; TDEF
ntim = 1260 ; time steps
nmos = 12
UNDEF = 1E+20
srcFileName = diri+fili
;---Output (and input) files
srcGridName = "ACCESS1-0_SCRIP.nc" ;
dstGridName = "World1deg_SCRIP.nc"
wgtFile = "ACCESS1-0_2_World_SCRIP.nc"
;---Retrieve data from BCC Grid
sfile = addfile(srcFileName,"r")
nlat = 180 ; YDEF
mlon = 360 ; XDEF
; tas_new = new((/ntim,nlat,mlon/),float,UNDEF)
; tos =new((/ntim,nlat,mlon/),float,UNDEF)
tas = sfile->tas(612:1871,:,:) ;
;---Set to True if you want to skip any of these steps
SKIP_BCC_SCRIP_GEN = False
SKIP_WORLD_SCRIP_GEN = False
SKIP_WGT_GEN = False
;----------------------------------------------------------------------
; Convert BCC to SCRIP file.
;----------------------------------------------------------------------
nav_lat = sfile->lat
nav_lon = sfile->lon
if(.not.SKIP_BCC_SCRIP_GEN) then
;---Convert to an SCRIP Convention file.
Opt = True
Opt@ForceOverwrite = True
Opt@PrintTimings = True
Opt@Mask2D = where(.not.ismissing(tas(0,:,:)),1,0)
rectilinear_to_SCRIP(srcGridName,nav_lat,nav_lon,Opt)
;---Clean up
delete(Opt)
end if
;----------------------------------------------------------------------
; Convert 1 degree world grid to SCRIP file
;----------------------------------------------------------------------
if(.not.SKIP_WORLD_SCRIP_GEN)
Opt = True
Opt@LLCorner = (/-89.75d, 0.00d /)
Opt@URCorner = (/ 89.75d, 359.75d /)
Opt@ForceOverwrite = True
Opt@PrintTimings = True
Opt@Title = "World grid 1x1 degree resolution"
latlon_to_SCRIP(dstGridName,"1deg",Opt)
;---Clean up
delete(Opt)
end if
;----------------------------------------------------------------------
; Generate interpolation weights for BCC Grid to World Grid
;----------------------------------------------------------------------
if(.not.SKIP_WGT_GEN) then
Opt = True
Opt@SrcESMF = False
Opt@DstESMF = False
Opt@ForceOverwrite = True
Opt@PrintTimings = True
ESMF_regrid_gen_weights(srcGridName, dstGridName, wgtFile, Opt)
;---Clean up
delete(Opt)
end if
;----------------------------------------------------------------------
; Interpolate data from BCC to World 1-degree grid.
;----------------------------------------------------------------------
Opt = True
; Opt@Debug = True
Opt@PrintTimings = True
tas_regrid = ESMF_regrid_with_weights(tas,wgtFile,Opt)
printVarSummary(tas_regrid)
tas_new = tas_regrid
printVarSummary(tas_new)
;---Clean up
delete(Opt)
;------output regridding data-------
NC = True
year = 1901 ; TDEF
ntim = 1260 ; time steps
nmos = 12
UNDEF = 1E+20
nlat = 180 ; YDEF
mlon = 360 ; XDEF
tas_new = new((/ntim,nlat,mlon/),float,UNDEF)
;=============================
; create coordinate variables
;============================
lon=fspan(0.75,359.75,360)
lon!0 = "lon"
lon@long_name = "longitude"
lon@units = "degrees_east"
lon&lon = lon
lat=fspan(-89.75,89.75,180)
lat!0 = "lat"
lat@long_name = "latitude"
lat@units = "degrees_north"
lat&lat = lat
years = sprinti("%4i",ispan(1901,2005,1)) ; write integers
months = sprinti("%0.2i", ispan(1,12,1)) ; ditto for
time = new(dimsizes(years)*dimsizes(months),"integer") ; create
k=0
do i=0,dimsizes(years)-1
do j = 0, dimsizes(months)-1
time(k)=stringtointeger(years(i)+months(j)) ; combine to
k=k+1
end do
end do
time!0 = "time"
time@long_name = "time"
time@units = "yyyymm"
time&time = time
tas_new!0 = "time"
tas_new!1 = "lat"
tas_new!2 = "lon"
tas_new&time = time
tas_new&lat = lat
tas_new&lon = lon
tas_new@long_name = "Near-Surface Temperature" ; Vars
tas_new@units = "K"
tas_new@_FillValue = 1E+20 ; UNDEF
printVarSummary(tas_new)
print ("min(tas_new)="+min(tas_new))
print ("max(tosregrid)="+max(tas_new))
TAS = new((/ntim,nlat,mlon/),float,UNDEF)
TAS(:,:,:) = tas_new(:,:,:)
copy_VarMeta(tas_new,TAS)
if (NC) then
diro = "./" ; Output directory
filo = "tas_Amon_ACCESS1-0_historical_r1i1p1_185001-201212-1x1.nc"
system("/bin/rm -f " + diro + filo) ; remove if exists
ncdf = addfile (diro+filo, "c")
; make time and UNLIMITED dimension ; recommended for most applications
filedimdef(ncdf,"time",-1,True)
; output variables directly
ncdf->TAS = TAS
end if
exit
end
(0) =====> CPU Elapsed Time: rectilinear_to_SCRIP: 3.63145 seconds <=====
(0) =====> CPU Elapsed Time: latlon_to_SCRIP: 8.18476 seconds <=====
(0) =====> CPU Elapsed Time: ESMF_regrid_gen_weights: 0.299954 seconds <=====
(0) ESMF_regrid_with_weights: warning: destination grid is not completely
(0) covered by the source grid.
(0) =====> CPU Elapsed Time: ESMF_regrid_with_weights: 14.9517 seconds <=====
The result is that there is no value in the new ......1x1.nc file except lat/lon information. I don't know why that warning: destination grid is not completely covered by the source grid. Actually, I don't know how to change Opt@Mask2D and Opt@LLCorner Opt@URCorner .
I change Opt@LLCorner to (/-90d, 0.00d /) and Opt@URCorner to (/ 90d, 360d /) but it can't work. Why? Can I use the Opt@Mask2D to mask ocean data or else?
NO.2(loop)
begin
;---Input file
diri = "./"
fili = "tas_Amon_ACCESS1-0_historical_r1i1p1_185001-200512.nc"
year = 1901 ; TDEF
ntim = 1260 ; time steps
nmos = 12
UNDEF = 1E+20
srcFileName = diri+fili
;---Output (and input) files
srcGridName = "ACCESS1-0_SCRIP.nc" ;
dstGridName = "World1deg_SCRIP.nc"
wgtFile = "ACCESS1-0_2_World_SCRIP.nc"
;---Retrieve data from BCC Grid
sfile = addfile(srcFileName,"r")
nlat = 180 ; YDEF
mlon = 360 ; XDEF
tas_new = new((/ntim,nlat,mlon/),float,UNDEF)
do nt = 0,ntim-1
tos = sfile->tas(nt+612,:,:) ;624 change to 612
;---Set to True if you want to skip any of these steps
SKIP_BCC_SCRIP_GEN = False
SKIP_WORLD_SCRIP_GEN = False
SKIP_WGT_GEN = False
;----------------------------------------------------------------------
; Convert BCC to SCRIP file.
;----------------------------------------------------------------------
nav_lat = sfile->lat
nav_lon = sfile->lon
if(.not.SKIP_BCC_SCRIP_GEN) then
;---Convert to an SCRIP Convention file.
Opt = True
Opt@ForceOverwrite = True
Opt@PrintTimings = True
Opt@Mask2D = where(.not.ismissing(tos),1,0)
rectilinear_to_SCRIP(srcGridName,nav_lat,nav_lon,Opt)
;---Clean up
delete(Opt)
end if
;----------------------------------------------------------------------
; Convert 1 degree world grid to SCRIP file
;----------------------------------------------------------------------
if(.not.SKIP_WORLD_SCRIP_GEN)
Opt = True
Opt@LLCorner = (/-89.75d, 0.00d /)
Opt@URCorner = (/ 89.75d, 359.75d /)
Opt@ForceOverwrite = True
Opt@PrintTimings = True
Opt@Title = "World grid 1x1 degree resolution"
latlon_to_SCRIP(dstGridName,"1deg",Opt)
;---Clean up
delete(Opt)
end if
;----------------------------------------------------------------------
; Generate interpolation weights for BCC Grid to World Grid
;----------------------------------------------------------------------
if(.not.SKIP_WGT_GEN) then
Opt = True
Opt@SrcESMF = False
Opt@DstESMF = False
Opt@ForceOverwrite = True
Opt@PrintTimings = True
ESMF_regrid_gen_weights(srcGridName, dstGridName, wgtFile, Opt)
;---Clean up
delete(Opt)
end if
;----------------------------------------------------------------------
; Interpolate data from BCC to World 1-degree grid.
;----------------------------------------------------------------------
Opt = True
; Opt@Debug = True
Opt@PrintTimings = True
tos_regrid = ESMF_regrid_with_weights(tos,wgtFile,Opt)
printVarSummary(tos_regrid)
tas_new(nt,:,:) = tos_regrid(:,:)
printVarSummary(tas_new)
;---Clean up
delete(Opt)
end do
;------output regridding data-------
NC = True
;=============================
; create coordinate variables
;============================
lon=fspan(0.75,359.75,360)
lon!0 = "lon"
lon@long_name = "longitude"
lon@units = "degrees_east"
lon&lon = lon
lat=fspan(-89.75,89.75,180)
lat!0 = "lat"
lat@long_name = "latitude"
lat@units = "degrees_north"
lat&lat = lat
years = sprinti("%4i",ispan(1901,2005,1)) ; write integers
months = sprinti("%0.2i", ispan(1,12,1)) ; ditto for
time = new(dimsizes(years)*dimsizes(months),"integer") ; create
k=0
do i=0,dimsizes(years)-1
do j = 0, dimsizes(months)-1
time(k)=stringtointeger(years(i)+months(j)) ; combine to
k=k+1
end do
end do
time!0 = "time"
time@long_name = "time"
time@units = "yyyymm"
time&time = time
tas_new!0 = "time"
tas_new!1 = "lat"
tas_new!2 = "lon"
tas_new&time = time
tas_new&lat = lat
tas_new&lon = lon
tas_new@long_name = "Near-Surface Temperature" ; Vars
tas_new@units = "K"
tas_new@_FillValue = 1E+20 ; UNDEF
printVarSummary(tas_new)
print ("min(tas_new)="+min(tas_new))
print ("max(tosregrid)="+max(tas_new))
TAS = new((/ntim,nlat,mlon/),float,UNDEF)
TAS(:,:,:) = tas_new(:,:,:)
copy_VarMeta(tas_new,TAS)
if (NC) then
diro = "./" ; Output directory
filo = "tas_Amon_ACCESS1-0_historical_r1i1p1_185001-201212-1x1.nc"
system("/bin/rm -f " + diro + filo) ; remove if exists
ncdf = addfile (diro+filo, "c")
; make time and UNLIMITED dimension ; recommended for most applications
filedimdef(ncdf,"time",-1,True)
; output variables directly
ncdf->TAS = TAS
end if
exit
Result: It seems to be correct but takes a long time (more than 3 hours).
> Subject: Re: About the new ESMF regridding !
> From: haley@ucar.edu
> Date: Fri, 2 Nov 2012 09:40:51 -0600
> CC: ncl-talk@ucar.edu
> To: chenxiaochen_ams@hotmail.com
>
> Dear Echo,
>
> Welcome to the NCL community.
>
> In the ESMF_all_7.ncl script, the "test_CNRM.nc" file actually only has one time step, so that's partially why it's only reading the one time step.
>
> However, if this file had multiple time steps, then you can easily change to multiple time steps by changing this line:
>
> tos = sfile->tos(0,:,:)
>
> to:
>
> tos = sfile->tos
>
>
> Assuming the mask is the same for each timestep, then you also need to change this line:
>
> Opt@Mask2D = where(.not.ismissing(tos),1,0)
>
> to:
>
> Opt@Mask2D = where(.not.ismissing(tos(0,:,:)),1,0)
>
> When you get to this call:
>
> tos_regrid = ESMF_regrid_with_weights(tos,wgtFile,Opt)
>
> It should automatically regrid across however many dimensions you have. "tos_regrid" will have the same leftmost dimensions as "tos", and the rightmost dimensions will be whatever grid you regridded it to.
>
> The important thing to note is that this method will only work if your grids are the same across all leftmost dimensions. If you have a source or destination grid that changes with each time or level or any other leftmost dimension, then you need to generate a separate weights file for each of these.
>
> If you look at example ESMF_all_6.ncl, which regrids data from a CMIP5 grid to a 1x1 grid, the same rules apply. This file has multiple timesteps and levels, so you would simply change two lines:
>
> thetao = sfile->thetao(0,0,:,:)
>
> to:
>
> thetao = sfile->thetao
>
> and:
>
> Opt@Mask2D = where(.not.ismissing(thetao),1,0)
>
> to:
>
> Opt@Mask2D = where(.not.ismissing(thetao(0,0,:,:)),1,0)
>
> In this case, then, since "thetao" is dimensioned 600 x 19 x 220 x 256 (time x lev x lat x lon), then "thetao_regrid" will be dimensioned 600 x 19 x 180 x 360, since 180 x 360 is the size of the destination (1x1) grid.
>
> Hope this helps clear things up.
>
> --Mary
>
>
> On Nov 1, 2012, at 8:03 PM, ³ÂÏþ³¿ wrote:
>
> > Dear ncl communit,
> >
> >
> > I am a new leaner. It's first time for me to send a email to the community. English is not my mother tougue, so maybe there are some mistakes in my email. Thank you for your understanding.
> >
> > I run the new ESMF regridding script- ESMF_all_7.ncl , and I succeed in plotting after I modify some instrctions.( http://www.ncl.ucar.edu/Applications/ESMF.shtml) I also scan other scripts in the ESMF regridding example pages for reference. But I want to regrid a 3D or 4D (or higher) data like A(:,:,:) and output it ,but all scripts just regrid a 2D data like A(0,:,:) or A(0,0,:,:).
> >
> > I have a large number of datum to regrid. My datum on rectilinear are from CMIP. The type of interpolation methods is bilinear.Do you have some suggestions for me? Maybe I should use cycling but I do not know how.
> >
> > Thanks for your attention!
> >
> >
> >
> > Echo
> > _______________________________________________
> > ncl-talk mailing list
> > List instructions, subscriber options, unsubscribe:
> > http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
Received on Sat Nov 3 03:55:02 2012
This archive was generated by hypermail 2.1.8 : Tue Nov 06 2012 - 15:05:49 MST