!=============================================================================== ! CVS $Id$ ! CVS $Source$ ! CVS $Name$ !=============================================================================== PROGRAM main use calendar_mod IMPLICIT NONE !--- includes --- #include !--- domain data --- integer :: ni ! size of i-axis of 2d domain integer :: nj ! size of j-axis of 2d domain integer :: nt ! time dimension real,allocatable :: xc( :) ! x-coords of center real,allocatable :: yc( :) ! y-coords of center real,allocatable :: data(:,: ) ! temporary data real,allocatable :: SST(:,:,:) ! SST real,allocatable :: ifrac(:,:,:) ! ice fraction real,allocatable :: xc2( :) ! x-coords of center - padded real,allocatable :: yc2( :) ! y-coords of center - padded real,allocatable :: ifrac2(:,:,:) ! ice fraction - padded real,allocatable :: SST2(:,:,:) ! SST - padded integer :: date (12) ! integer date integer :: datesec(12) ! elapsed secs on date real :: time (12) ! elapsed days since ref date integer :: fdate (12) ! time coordinate (date) from data file character :: c !--- 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 :: yy,mm,dd ! year & month indicies integer :: i,j,k,n,m ! generic indicies integer :: beg(3) ! nc hyperslab begin vector integer :: cnt(3) ! nc hyperslab count vector character,allocatable :: strarr(:) ! variable length char string character(LEN=256) :: 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 character(LEN=2048) :: strlong !--- 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 read Hurrell/Hadley/STR SST+iceFraction data file and format appropriately ! for dice5 & docn5 ! ! NOTES: ! o to compile on an NCAR's SGI, utefe (Jul 2001): ! unix> f90 -64 -mips4 -O2 -r8 -i4 -I /usr/local/include \ ! -L/usr/local/lib64/r4i4 -lnetcdf hurrellData.F90 !------------------------------------------------------------------------------- !---------------------------------------------------------------------------- write(6,F00) 'create dice5/docn5 sst/ifrac data from original hurrell data' !---------------------------------------------------------------------------- fn_in = '/fs/cgd/data0/shea/jhurrell/MODEL.had+oiv2.sst.mnly.49-01.unf.nc' fn_out = 'sst_clim_hurrell_1x1_yymmdd.nc' write(*,*) 'fn_in = ',fn_in (1:len_trim(fn_in )) !---------------------------------------------------------------------------- write(6,F00) 'input SST data...' !---------------------------------------------------------------------------- !--- document global attributes ----------------- write(6,F00) 'o file = ',fn_in(1:len_trim(fn_in)) rcode = nf_open(trim(fn_in),NF_NOWRITE,fid) if (rcode /= NF_NOERR) write(*,F00) nf_strerror(rcode) !---------------------------------------------- ! get domain info for source grid !---------------------------------------------- rcode = nf_inq_dimid (fid, 'lon' ,did) rcode = nf_inq_dimlen(fid, did , ni) rcode = nf_inq_dimid (fid, 'lat' ,did) rcode = nf_inq_dimlen(fid, did , nj) rcode = nf_inq_dimid (fid, 'time',did) rcode = nf_inq_dimlen(fid, did , nt) write(6,F02) 'o ni,nj,nt=',ni,nj,nt allocate( xc (ni ) ) ! x-coordinates of center allocate( yc ( nj) ) ! y-coordinates of center allocate( data (ni,nj) ) ! temporary data input allocate( SST (ni,nj,12)) ! sst climatology allocate( ifrac(ni,nj,12)) ! ifrac climatology allocate( xc2 (0:ni+1 ) ) ! x-coordinates of center - padded allocate( yc2 (0:nj+1) ) ! y-coordinates of center - padded allocate( SST2 (0:ni+1,0:nj+1,12)) ! sst climatology - padded allocate( ifrac2(0:ni+1,0:nj+1,12)) ! ifrac climatology - padded rcode = nf_inq_varid (fid,'lon',vid) rcode = nf_get_var_double(fid,vid, xc ) rcode = nf_inq_varid (fid,'lat',vid) rcode = nf_get_var_double(fid,vid, yc ) !---------------------------------------------- ! time,date, & sec for monthly climatology !---------------------------------------------- rcode = nf_inq_varid (fid,'time',vid) rcode = nf_get_vara_int(fid,vid,1,12,fdate) do i=1,12 call date2ymd(fdate(i),yy,mm,dd) call ymd2eday( 0,mm,dd,n) ! n is elapsed day since start of year time (i) = n date (i) = mm*100 + dd datesec(i) = 0 end do !---------------------------------------------- ! make monthly SST climatology !---------------------------------------------- rcode = nf_inq_varid (fid,'SST' ,vid) beg(1) = 1 beg(2) = 1 beg(3) = 0 cnt(1) = ni cnt(2) = nj cnt(3) = 1 if (nt /= (2001-1949+1)*12) then write(*,*) 'ERROR: what dates are in data file?' stop end if SST = 0.0 n=0 do yy=1949,2001 !! do yy=1949,1951 do mm=1,12 if (mm==1) write(*,*) 'reading SST: year = ',yy n=n+1 beg(3) = n rcode = nf_get_vara_double(fid,vid,beg,cnt,data) SST(:,:,mm) = SST(:,:,mm) + data end do end do !! SST = SST/float(1951-1949+1) SST = SST/float(2001-1949+1) write(*,*) 'Min/Max Val SST = ',minval(SST),maxval(SST) rcode = nf_close(fid) !---------------------------------------------------------------------------- write(6,F00) 'input ifrac data...' !---------------------------------------------------------------------------- fn_in = '/fs/cgd/data0/shea/jhurrell/MODEL.had+oiv2.ice.mnly.49-01.unf.nc' write(6,F00) 'o file = ',trim(fn_in) rcode = nf_open(trim(fn_in),NF_NOWRITE,fid) if (rcode /= NF_NOERR) write(*,F00) nf_strerror(rcode) !---------------------------------------------- ! make monthly ifrac climatology !---------------------------------------------- rcode = nf_inq_varid(fid,'SEAICE' ,vid) beg(1) = 1 beg(2) = 1 beg(3) = 0 cnt(1) = ni cnt(2) = nj cnt(3) = 1 ifrac = 0.0 n=0 do yy=1949,2001 !! do yy=1949,1951 do mm=1,12 if (mm==1) write(*,*) 'reading ifrac: year = ',yy n=n+1 beg(3) = n rcode = nf_get_vara_double(fid,vid,beg,cnt,data) ifrac(:,:,mm) = ifrac(:,:,mm) + data end do end do !! ifrac = .01*ifrac/float(1951-1949+1) ifrac = .01*ifrac/float(2001-1949+1) ! convert from % to fraction write(*,*) 'Min/Max Val ifrac = ',minval(ifrac),maxval(ifrac) rcode = nf_close(fid) !---------------------------------------------------------------------------- write(6,F00) 'create new padded data' !---------------------------------------------------------------------------- xc2(0 ) = xc(ni) - 360.0 xc2(ni+1) = xc( 1) + 360.0 xc2(1:ni) = xc(1:ni) yc2(0 ) = -90.0 yc2(nj+1) = 90.0 yc2(1:nj) = yc(1:nj) do k = 1,12 !--- north/south pole padding --- SST2 (:, 0,k) = SUM(SST (:, 1,k))/float(ni) ! south-pole value SST2 (:,nj+1,k) = SUM(SST (:,nj,k))/float(ni) ! north-pole value ifrac2(:, 0,k) = SUM(ifrac(:, 1,k))/float(ni) ! south-pole value ifrac2(:,nj+1,k) = SUM(ifrac(:,nj,k))/float(ni) ! north-pole value !--- east/west periodic padding --- SST2 ( 0,1:nj,k) = SST (ni,:,k) ! wrap-around data SST2 (ni+1,1:nj,k) = SST (1 ,:,k) ! wrap-around data ifrac2( 0,1:nj,k) = ifrac(ni,:,k) ! wrap-around data ifrac2(ni+1,1:nj,k) = ifrac(1 ,:,k) ! wrap-around data end do SST2 (1:ni,1:nj,:) = SST (:,:,:) ifrac2(1:ni,1:nj,:) = ifrac(:,:,:) !---------------------------------------------------------------------------- write(6,F00) 'create new netCDF file...' !---------------------------------------------------------------------------- rcode = nf_create(trim(fn_out),NF_CLOBBER,fid) !---------------------------------------------------------------------------- !---------------------------------------------------------------------------- c = char(10) strlong = & & "Jim Hurrell (NCAR) combined the Hadley Center anomalies " // c // & & "with the Reynolds SST climatology to produce SSTs over " // c // & & "the ocean only (2001?). Several tests were run to determine " // c // & & "how values over land were to be created. " // c // & & " /fs/cgd/home0/shea/ncld/jhurrell/sst.ncl " // c // & & "Tests indicated that direct interpolation of the full SST " // c // & & "grids would be satisfactory. Note: values over land " // c // & & "are entirely bogus. In particular, values over Siberia " // c // & & "are particularly unusual. However, values near " // c // & & "coasts are close to those at nearby oceans. " // c // & & " " // c // & & "The grids produced for this file were created using " // c // & & "the NCL function cssgrid which interpolates randomly " // c // & & "spaced data on a sphere to a lat/lon grid. The original " // c // & & "grid values were then used to overwrite the interpolated " // c // & & "values over the oceans. The land values are the ONLY " // c // & & "interpolated values on these grids. The NCL script is: " // c // & & " /fs/cgd/home0/shea/ncld/jhurrell/sstDirect.ncl " // c // & & " " // c // & & "/fs/cgd/data0/shea/jhurrell/had+oiv2.sst.mnly.49-01.unf.nc " // c // & & "This nc file was created from the binary source file using: " // c // & & " /fs/cgd/home0/shea/ncld/jhurrell/sstBinary2nc.ncl " // c // & & " ------------------------- " // c // & & "Jim Hurrell provided the sea-ice percentages. " // c // & & "Basically, linear interpolation [East-west] was used. " // c // & & " /fs/cgd/home0/shea/ncld/jhurrell/iceDirect.ncl " // c // & & " " // c // & & "/fs/cgd/data0/shea/jhurrell/had+oiv2.ice.mnly.49-01.unf.nc " // c // & & "This nc file was created from the binary source file using: " // c // & & " /fs/cgd/home0/shea/ncld/jhurrell/iceBinary2nc.ncl " // c // & & " ------------------------- " // c // & & "Brian Kauffman (NCAR) formed a 12 month climatology and combined "// c // & & "the SST and ice fraction data (described above) into one file " // c // & & "appropriate for used by the docn5 & dice5 CCSM2 data models. "// c // & & "The data were padded with -0.5 & 360.5 longitude bands " // c // & & "plus -90 and +90 latitude bands to facilitate interpolation " // c // & & "onto other grids (eg. pop shifted pole, etc.). " !---------------------------------------------------------------------------- write(6,F00) 'add global attributes' !---------------------------------------------------------------------------- str = 'Hurrell SST climatology (Reynolds SST + Hadley anomolies)' strlen=len_trim(str) rcode = nf_put_att_text(fid,NF_GLOBAL,'title' ,strlen,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) strlen=len_trim(str) rcode = nf_put_att_text(fid,NF_GLOBAL,'history' ,strlen,str) str ='CF-1.0' strlen=len_trim(str) rcode = nf_put_att_text(fid,NF_GLOBAL,'Conventions',strlen,str) str='/fs/cgd/data0/shea/jhurrell/MODEL.had+oiv2.sst.mnly.49-01.unf.nc '//c//& & '/fs/cgd/data0/shea/jhurrell/MODEL.had+oiv2.ice.mnly.49-01.unf.nc ' strlen=len_trim(str) rcode = nf_put_att_text(fid,NF_GLOBAL,'source' ,strlen,str) strlen=len_trim(strlong) rcode = nf_put_att_text(fid,NF_GLOBAL,'comments' ,strlen,strlong) !---------------------------------------------------------------------------- write(6,F00) 'define dimension data' !---------------------------------------------------------------------------- rcode = nf_def_dim(fid, 'lon' , ni+2, did) rcode = nf_def_dim(fid, 'lat' , nj+2, did) rcode = nf_def_dim(fid, 'time', 12, did) rcode = nf_inq_dimid(fid,'lon' ,vdid(1)) rcode = nf_inq_dimid(fid,'lat' ,vdid(2)) rcode = nf_inq_dimid(fid,'time',vdid(3)) rcode = nf_def_var (fid,'lon',NF_FLOAT ,1,vdid(1),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,'lat',NF_FLOAT ,1,vdid(2),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_def_var(fid,'time',NF_REAL,1,did,vid) str = 'time' rcode = nf_put_att_text(fid,vid,'long_name',len_trim(str),str) str = 'days since 0000-01-01 00:00:00' rcode = nf_put_att_text(fid,vid,'units' ,len_trim(str),str) str = 'noleap' rcode = nf_put_att_text(fid,vid,'calendar' ,len_trim(str),str) !---------------------------------------------------------------------------- write(6,F00) 'define field data' !---------------------------------------------------------------------------- rcode = nf_def_var (fid,'date',NF_INT ,1,vdid(3),vid) str = 'calendar date' rcode = nf_put_att_text(fid,vid,"long_name",len_trim(str),str) str = 'yyyymmdd' rcode = nf_put_att_text(fid,vid,"units" ,len_trim(str),str) rcode = nf_def_var (fid,'datesec',NF_INT ,1,vdid(3),vid) str = 'seconds elapsed on calendar date' rcode = nf_put_att_text(fid,vid,"long_name",len_trim(str),str) str = 'seconds' rcode = nf_put_att_text(fid,vid,"units" ,len_trim(str),str) rcode = nf_def_var (fid,'sst',NF_FLOAT ,3,vdid,vid) str = 'sea surface temperature' rcode = nf_put_att_text(fid,vid,"long_name",len_trim(str),str) str = 'degrees C' rcode = nf_put_att_text(fid,vid,"units" ,len_trim(str),str) rcode = nf_def_var (fid,'ifrac',NF_FLOAT ,3,vdid,vid) str = 'ice fraction' rcode = nf_put_att_text(fid,vid,"long_name",len_trim(str),str) str = 'fraction ranges from 0 to 1' rcode = nf_put_att_text(fid,vid,"comment" ,len_trim(str),str) rcode = nf_enddef(fid) !---------------------------------------------------------------------------- write(6,F00) 'put field data' !---------------------------------------------------------------------------- rcode = nf_inq_varid (fid, 'lon' , vid) rcode = nf_put_var_double(fid, vid , xc2) rcode = nf_inq_varid (fid, 'lat' , vid) rcode = nf_put_var_double(fid, vid , yc2) rcode = nf_inq_varid (fid, 'sst' , vid) rcode = nf_put_var_double(fid, vid , sst2) rcode = nf_inq_varid (fid, 'ifrac' , vid) rcode = nf_put_var_double(fid, vid , ifrac2) rcode = nf_inq_varid (fid, 'time' , vid) rcode = nf_put_var_double(fid, vid , time) rcode = nf_inq_varid (fid, 'date' , vid) rcode = nf_put_var_int (fid, vid , date) rcode = nf_inq_varid (fid, 'datesec', vid) rcode = nf_put_var_int (fid, vid , datesec) !---------------------------------------------------------------------------- write(6,F00) 'Close file...' !---------------------------------------------------------------------------- rcode = nf_close(fid) END PROGRAM !===============================================================================