!=============================================================================== ! CVS $Id$ ! CVS $Source$ ! CVS $Name$ !=============================================================================== PROGRAM make_domain IMPLICIT NONE !--- includes --- include 'netcdf.inc' ! netCDF defs !--- domain data --- integer :: n ! size of 1d domain integer :: ni ! size of i-axis of 2d domain integer :: nj ! size of j-axis of 2d domain integer :: nv = 4 ! assume retalinear grid real ,pointer :: xc( :) ! x-coords of center real ,pointer :: yc( :) ! y-coords of center real ,pointer :: xv(:,:) ! x-coords of verticies real ,pointer :: yv(:,:) ! y-coords of verticies integer,pointer :: mask( :) ! domain mask real ,pointer :: area( :) ! cell area !--- local --- character(LEN= 8) :: cdate ! wall clock date character(LEN=10) :: ctime ! wall clock time integer :: rcode ! routine return error code integer :: fid ! nc file ID integer :: vid ! nc variable ID integer :: did ! nc dimension ID integer :: vdid(3) ! vector of nc dimension ID integer :: i,j,k,n,m ! generic indicies character,allocatable :: strarr(:) ! variable length char string character(LEN=80) :: str ! fixed length char string character(LEN=80) :: str_title ! global attribute str - title character(LEN=80) :: str_source ! global attribute str - source character(LEN=80) :: source ! optional global att str - source character(LEN=80) :: title ! replacement global att str - title integer :: strlen ! (trimmed) length of string character(LEN=80) :: fn_in ! file name ( input nc file) character(LEN=80) :: fn_out ! file name (output nc file) character(LEN=32) :: attstr ! netCDF attribute name string real,parameter :: pi = 3.14159265358979323846 NAMELIST / in_parm / fn_in, fn_out, title, source !--- formats --- character(LEN=*),parameter :: F00 = "(120a )" character(LEN=*),parameter :: F02 = "(a,4i6)" character(LEN=*),parameter :: F10=& & "('Data created: 'i4,'-',i2,2('-',i2),' ',i2,2('-',i2),' ')" !------------------------------------------------------------------------------- ! PURPOSE: ! o given a SCRIP map matrix data file, create a docn/dice domain data file ! ! NOTES: ! o all output data is base on the "_a" grid, the "_b" grid is ignored ! o to compile on an NCAR's SGI, utefe (Feb 2001): ! unix> f90 -64 -mips4 -r8 -i4 -lfpe -I/usr/local/include Make_domain.F90 \ ! -L/usr/local/lib64/r4i4 -lnetcdf !------------------------------------------------------------------------------- write(6,F00) 'create a dice5/docn5 domain file from a scrip matrix data file' !---------------------------------------------------------------------------- write(6,F00) 'input namelist data...' !---------------------------------------------------------------------------- fn_in ='null' fn_out='null' title ='null' source='null' read(*,nml=in_parm) write(*,*) 'fn_in = ',fn_in (1:len_trim(fn_in )) write(*,*) 'fn_out = ',fn_out(1:len_trim(fn_out)) write(*,*) 'title = ',title (1:len_trim(title )) write(*,*) 'source = ',source(1:len_trim(source)) !---------------------------------------------------------------------------- write(6,F00) 'input SCRIP data...' !---------------------------------------------------------------------------- !--- document global attributes ----------------- write(6,F00) 'o file = ',fn_in(1:len_trim(fn_in)) rcode = nf_open(fn_in(1:len_trim(fn_in)),NF_NOWRITE,fid) if (rcode.ne.NF_NOERR) write(*,F00) nf_strerror(rcode) rcode = nf_get_att_text(fid, NF_GLOBAL, 'title' , str_title ) write(6,F00) 'o title = ',str_title (1:len_trim(str_title )) if (title(1:4) /= 'null') str_title=title rcode = nf_get_att_text(fid, NF_GLOBAL, 'source' , str_source) write(6,F00) 'o source = ',str_source(1:len_trim(str_source)) if (source(1:4) /= 'null') str_source=source !---------------------------------------------- ! get domain info for source grid !---------------------------------------------- rcode = nf_inq_dimid (fid, 'n_a' ,did) rcode = nf_inq_dimlen(fid, did , n) rcode = nf_inq_dimid (fid, 'nv_a',did) rcode = nf_inq_dimlen(fid, did , nv) rcode = nf_inq_dimid (fid, 'ni_a',did) rcode = nf_inq_dimlen(fid, did , ni) rcode = nf_inq_dimid (fid, 'nj_a',did) rcode = nf_inq_dimlen(fid, did , nj) write(6,F02) 'o n,ni,nj,nv=',n,ni,nj,nv allocate( xc( n)) ! x-coordinates of center allocate( yc( n)) ! y-coordinates of center allocate( xv(nv,n)) ! x-coordinates of verticies allocate( yv(nv,n)) ! y-coordinates of verticies allocate(mask( n)) ! domain mask allocate(area( n)) ! grid cell area rcode = nf_inq_varid (fid,'xc_a',vid) rcode = nf_get_var_double(fid,vid, xc ) rcode = nf_inq_varid (fid,'yc_a',vid) rcode = nf_get_var_double(fid,vid, yc ) rcode = nf_inq_varid (fid,'xv_a',vid) rcode = nf_get_var_double(fid,vid, xv ) rcode = nf_inq_varid (fid,'yv_a',vid) rcode = nf_get_var_double(fid,vid, yv ) rcode = nf_inq_varid (fid,'mask_a',vid ) rcode = nf_get_var_int (fid,vid,mask ) rcode = nf_inq_varid (fid,'area_a',vid ) rcode = nf_get_var_double(fid,vid,area ) rcode = nf_close(fid) ! !---------------------------------------------------------------------------- write(6,F00) 'output CSM data...' !---------------------------------------------------------------------------- !----------------------------------------------------------------- ! create a new nc file !----------------------------------------------------------------- rcode = nf_create(fn_out(1:len_trim(fn_out)),NF_CLOBBER,fid) if (rcode.ne.NF_NOERR) write(*,F00) nf_strerror(rcode) !----------------------------------------------------------------- ! global attributes !----------------------------------------------------------------- str = 'domain data: '//str_title rcode = nf_put_att_text(fid,NF_GLOBAL,'title' ,len_trim(str),str) str = 'CCSM2.0 data model domain description' rcode = nf_put_att_text(fid,NF_GLOBAL,'conventions',len_trim(str),str) str = source if ( str(1:4) /= 'null' ) & & rcode = nf_put_att_text(fid,NF_GLOBAL,'source' ,len_trim(str),str) call date_and_time(cdate,ctime) str = 'File created: '//cdate(1:4)//'-'//cdate(5:6)//'-'//cdate(7:8) & & //' '//ctime(1:2)//':'//ctime(3:4)//':'//ctime(5:6) rcode = nf_put_att_text(fid,NF_GLOBAL,'history' ,len_trim(str),str) !----------------------------------------------------------------- ! dimension data !----------------------------------------------------------------- if (n /= ni*nj) STOP 'n' rcode = nf_def_dim(fid, 'n' , n , did) ! # of points total rcode = nf_def_dim(fid, 'ni', ni, did) ! # of points wrt i rcode = nf_def_dim(fid, 'nj', nj, did) ! # of points wrt j rcode = nf_def_dim(fid, 'nv', 4, did) ! # of verticies per cell !----------------------------------------------------------------- ! define data -- coordinates, input grid !----------------------------------------------------------------- rcode = nf_inq_dimid(fid,'n' , did ) rcode = nf_inq_dimid(fid,'ni',vdid(1)) rcode = nf_inq_dimid(fid,'nj',vdid(2)) rcode = nf_def_var (fid,'xc',NF_DOUBLE,2,vdid,vid) str = 'longitude of grid cell center' rcode = nf_put_att_text(fid,vid,"long_name",len_trim(str),str) str = 'degrees east' rcode = nf_put_att_text(fid,vid,"units" ,len_trim(str),str) rcode = nf_def_var (fid,'yc',NF_DOUBLE,2,vdid,vid) str = 'latitude of grid cell center' rcode = nf_put_att_text(fid,vid,"long_name",len_trim(str),str) str = 'degrees north' rcode = nf_put_att_text(fid,vid,"units" ,len_trim(str),str) rcode = nf_inq_dimid(fid,'nv',vdid(1)) rcode = nf_inq_dimid(fid,'ni',vdid(2)) rcode = nf_inq_dimid(fid,'nj',vdid(3)) rcode = nf_def_var (fid,'xv',NF_DOUBLE,3,vdid,vid) str = 'longitude of grid cell verticies' rcode = nf_put_att_text(fid,vid,"long_name",len_trim(str),str) str = 'degrees east' rcode = nf_put_att_text(fid,vid,"units" ,len_trim(str),str) rcode = nf_def_var (fid,'yv',NF_DOUBLE,3,vdid,vid) str = 'latitude of grid cell verticies' rcode = nf_put_att_text(fid,vid,"long_name",len_trim(str),str) str = 'degrees north' rcode = nf_put_att_text(fid,vid,"units" ,len_trim(str),str) rcode = nf_inq_dimid(fid,'ni',vdid(1)) rcode = nf_inq_dimid(fid,'nj',vdid(2)) rcode = nf_def_var (fid,'mask',NF_INT ,2,vdid,vid) str = 'domain mask' rcode = nf_put_att_text(fid,vid,"long_name",len_trim(str),str) str = 'unitless' rcode = nf_put_att_text(fid,vid,"units" ,len_trim(str),str) rcode = nf_def_var (fid,'area',NF_DOUBLE,2,vdid,vid) str = 'area of grid cell in radians squared' rcode = nf_put_att_text(fid,vid,"long_name",len_trim(str),str) str = 'area' rcode = nf_put_att_text(fid,vid,"units" ,len_trim(str),str) !----------------------------------------------------------------- write(6,F00) 'o put data...' !----------------------------------------------------------------- rcode = nf_enddef(fid) rcode = nf_inq_varid (fid, 'xc',vid) rcode = nf_put_var_double(fid, vid , xc) rcode = nf_inq_varid (fid, 'yc',vid) rcode = nf_put_var_double(fid, vid , yc) rcode = nf_inq_varid (fid, 'xv',vid) rcode = nf_put_var_double(fid, vid , xv) rcode = nf_inq_varid (fid, 'yv',vid) rcode = nf_put_var_double(fid, vid , yv) rcode = nf_inq_varid (fid,'mask',vid) rcode = nf_put_var_int (fid, vid ,mask) rcode = nf_inq_varid (fid,'area',vid) rcode = nf_put_var_double(fid, vid ,area) rcode = nf_close(fid) if (rcode.ne.NF_NOERR) write(*,F00) nf_strerror(rcode) END PROGRAM !===============================================================================