PROGRAM NETCDFFROMWOABIN_OHCNC4 c@***h* OCEAN_HEAT_CONTENT/netcdffromwoabin_ohcnc4.f c c NAME c netcdffromwoabin_ohcnc4.f c c PURPOSE c Specifically for OCEAN_HEAT_CONTENT (OHC), the gridded c temperature anomaly fields and related statistics, and c the ocean heat content fields and integrated time series c are placed in a Climate-Forecast (CF) convention netCDF c file c DESCRIPTION c netcdffromwoabin_ohcnc4.f is a multipurpose program within c the World Ocean Database (WOD) system. It reads in native c WOD binary format and writes out in a CF compliant netCDF c file. c c INPUT c file woatoncfile.d lists the files to convert from WOD c binary form to CF compliant netCDF and the the name of c the variable, units, grid size, temporal span (e.g. seasonal), c time period (e.g. 1955 to 2016), geographical area (global c for OHC), and grid size (1-degree for OHC). woatoncfile.d c also lists all integrated time series data files to also c include. File global_attributes.anomaly.d provides input c information for the global attributes field of the netCDF c file. File woa_field_definitionsnc.txt provides full definitions c of all the statistical fields which can be included in c the netCDF file, which are requested in woatoncfile.d by c code number (first column in woa_field_definitionsnc.dat). c File timeseries_area_info.txt provides naming conventions for c ocean basins, also requested in woatoncfile.d by the code in c the first column. The input data files, as given in c woatoncfile.d are all in WOD binary form, excepting the c time series integrated files which are in ASCII. c c OUTPUT c The output files, whose name and directory are given in woatoncfile.d c are CF compliant netCDF files of all fields designated in c woatoncfile.d, self-described. c cAUTHOR c Tim Boyer, National Centers for Environmental Information (NCEI) c cCOPYRIGHT cTHIS SOFTWARE AND ITS DOCUMENTATION ARE CONSIDERED TO BE IN THE cPUBLIC DOMAIN AND THUS ARE AVAILABLE FOR UNRESTRICTED PUBLIC USE. cTHEY ARE FURNISHED "AS IS." THE AUTHORS, THE UNITED STATES GOVERNMENT, cITS INSTRUMENTALITIES, OFFICERS, EMPLOYEES, AND AGENTS MAKE NO cWARRANTY, EXPRESS OR IMPLIED, AS TO THE USEFULNESS OF THE SOFTWARE cAND DOCUMENTATION FOR ANY PURPOSE. THEY ASSUME NO RESPONSIBILITY c(1) FOR THE USE OF THE SOFTWARE AND DOCUMENTATION; OR (2) TO PROVIDE cTECHNICAL SUPPORT TO USERS. c cHISTORY c August 14, 2017 Initial version prepared for Reference Environmental c Data Record (REDR) program c c***** C NETCDFFROMWOABIN CREATES A NETCDF FILE FROM WOA C BINARY FORM parameter (idimax=3600,jdimax=1800) parameter (kdimax=137, maxvars=100) parameter (intfill=-32767,imiss=-100) parameter (rfill=9.96921e+36,bmiss=-1.E10) parameter (maxts=50,maxyear=12000) character*80 varname,varnamefull,cunit character*80 varnamets,varnamefullts,cunitts character*200 infilename(maxvars),outfilename,landmask character*200 filename,infilenamets(maxts) character*80 globalfile character*500 ftitle dimension x(idimax,jdimax),dz(kdimax) dimension ival(idimax,jdimax),m1(idimax,jdimax) dimension ivars(maxvars),ivartsid(maxts) dimension ivardims(4),ivardimsize(4),ivarid(maxvars) dimension ivarindex(4),idz(kdimax) dimension intfile(maxvars) dimension nfieldtypets(maxts),icolumnts(maxts) dimension timeout(maxyear),iserrtsid(maxts) include 'netcdf.inc' infile=-1 c Get information for file call getinputinfnc(idim,jdim,xgrid,levels,dz, * kdimax,idz,nvars,maxvars,ivars,varname,istart,iend, * jstart,jend,ntime,ntimestart,ntimestart2,ntimeinterval, * infilename, outfilename,ncvar,cunit,nunit,globalfile, * ftitle,ifsize,nfirstyear,ntimeend,ntimeend2, * varnamefull,ncvarf,ntseries,nfieldtypets, * varnamets,nspvts,cunitts,nsputs,infilenamets, * varnamefullts,nspvfts,maxts,icolumnts,landmask,noland) c See if all fields are integrals iallintegral=0 do 25 nv=1,nvars call setupfieldvarnc(ivars(nv),varname,ncvar,cunit,nunit, * ivarid(nv),ivardims,ivardimsize,ncid,intfill,rfill, * intfile(nv),varnamefull,ncvarf,iallintegral,1) 25 continue if ( iallintegral .eq. nvars ) then iallintegral=1 levels0=1 else iallintegral=0 levels0=levels endif c open and read land mask if necessary landfile=-1 if ( noland .eq. 0 ) then call addCend(landmask,200,0) call fileopen(landmask,'rb+'//CHAR(0),landfile) call readleveli(1,1,1,1,1,idim,jdim,idim,jdim,idimax,jdimax, * 1,1,landfile,m1) call onefileclose(landfile) endif c set up netcdf file call fileopennc_classic(ncid,outfilename) call setupcrsnc(ncid,icrsid,0) instat=nf_def_dim(ncid,'nbounds',2,nbounddimid) if ( instat .ne. 0 ) then write(6,*) 'Problem assigning bounds dimension',instat endif jspan=(jend-jstart+1) instat=nf_def_dim(ncid,'lat',jspan,jdimid) if ( instat .ne. 0 ) then write(6,*) 'Problem assigning latitude dimension',instat endif ivardimsize(2)=jspan call latitudesetup(ncid,jdimid,ilatid,nbounddimid, * ilatboundid) if ( iend .gt. istart ) then ispan=iend-istart+1 else ispan=((idim-istart)+iend)+1 endif c write(6,*) 'ispan',idim,istart,iend,ispan instat=nf_def_dim(ncid,'lon',ispan,idimid) if ( instat .ne. 0 ) then write(6,*) 'Problem assigning longitude dimension',instat endif ivardimsize(1)=ispan call longitudesetup(ncid,idimid,ilonid,nbounddimid, * ilonboundid) c write(6,*) 'lonbound',ilonboundid instat=nf_def_dim(ncid,'depth',levels0,izdimid) if ( instat .ne. 0 ) then write(6,*) 'Problem assigning depth dimension',instat endif ivardimsize(3)=levels0 call depthsetup(ncid,izdimid,izid,nbounddimid,izboundid, * iallintegral) instat=nf_def_dim(ncid,'time',ntime,itimedimid) if ( instat .ne. 0 ) then write(6,*) 'Problem assigning time dimension',instat endif ivardimsize(4)=ntime call timesetup(ncid,itimedimid,itimeid,nbounddimid, * itimeboundid,nfirstyear) ivardims(4)=itimedimid ivardims(3)=izdimid ivardims(2)=jdimid ivardims(1)=idimid write(6,*) 'vardims',ivardims write(6,*) 'vardimsize',ivardimsize do 30 nv=1,nvars call setupfieldvarnc(ivars(nv),varname,ncvar,cunit,nunit, * ivarid(nv),ivardims,ivardimsize,ncid,intfill,rfill, * intfile(nv),varnamefull,ncvarf,iallintegral,0) 30 continue c write(6,*) 'setup',ntseries,iallintegral,nvars,itimeid call setuptimeseriesnc(ntseries,nfieldtypets, * varnamefullts,nspvfts,cunitts,nsputs,varnamets,nspvts, * ivartsid,iserrtsid,maxts,itimedimid,ncid,idimid,jdimid, * ibasinid) c Set up global variables c write(6,*) 'setglob' call setup_globalatt_woanc(ncid,globalfile,jstart,jend,istart, * iend,xgrid,dz(1),dz(levels),ftitle,ifsize,ntimestart, * ntimestart2,ntimeinterval,ntimeend,ntimeend2,outfilename) c write(6,*) 'setglob' instat=nf_enddef(ncid) if ( instat .ne. 0 ) then write(6,*) 'enddef not executed properly',instat endif c Enter latitudes, bounds call putlatitudenc(ncid,jstart,jend,xgrid,ilatid, * ilatboundid) c Enter longitudes, bounds c write(6,*) 'lonbound2',ilonboundid call putlongitudenc(ncid,istart,iend,xgrid,ilonid, * ilonboundid) c Enter time periods, bounds call puttimenc(ncid,ntime,ntimestart,ntimestart2, * ntimeinterval,itimeid,itimeboundid,nfirstyear, * timeout,maxyear) c Enter depths, bounds call putdepthnc(ncid,levels,dz,kdimax,izid, * izboundid,iallintegral,levels0) do 50 itz=1,ntime ivarindex(4)=itz do 55 nv=1,nvars call setwoafilename(filename,ntimestart,ntimeinterval, * ntimestart2,infilename(nv),itz) call fileopen(filename,'rb+'//CHAR(0),infile) do 57 kx=1,levels0 ivarindex(3)=idz(kx) if ( intfile(nv) .eq. 0 .or. ivars(nv) .eq. 7) then call readlevel(1,1,1,1,1,idim,jdim,idim,jdim,idimax,jdimax, * idz(kx),idz(kx),infile,x) else call readleveli(1,1,1,1,1,idim,jdim,idim,jdim,idimax,jdimax, * idz(kx),idz(kx),infile,ival) endif c write(6,*) 'kx',kx,x(181,151),ival(181,151) do 60 j=jstart,jend ivarindex(2)=(j-jstart+1) itimex=1 if ( istart .gt. iend ) itimex=2 index1=0 do 62 itx=1,itimex istart1=istart iend1=iend if ( itimex .gt. 1 ) then if ( itx .eq. 1 ) iend1=idim if ( itx .eq. 2 ) istart1=1 endif do 65 i=istart1,iend1 if ( noland .eq. 0 .and. m1(i,j) .le. idz(kx) ) then x(i,j)=bmiss ival(i,j)=imiss endif index1=index1+1 ivarindex(1)=index1 if ( intfile(nv) .eq. 0 ) then c write(6,*) 'nv',nv,ivarindex,i,j,x(i,j) if ( x(i,j) .le. bmiss ) x(i,j)=rfill instat=nf_put_var1_real(ncid,ivarid(nv),ivarindex,x(i,j)) else c Grid points data (7) if ( ivars(nv) .eq. 7 ) then if ( x(i,j) .gt. bmiss ) then ival(i,j)=x(i,j) else ival(i,j)=imiss endif endif if ( ival(i,j) .le. imiss ) ival(i,j)=intfill instat=nf_put_var1_int(ncid,ivarid(nv),ivarindex,ival(i,j)) endif if ( instat .ne. 0 ) then write(6,*) 'Problem writing out value',instat,x(i,j), * ival(i,j),nv,ivarindex,i,j endif 65 continue 62 continue 60 continue 57 continue call onefilecloseft(infile) 55 continue 50 continue c Basin mask if ( ntseries .gt. 0 ) then call putbasindefs(ncid,ibasinid,jstart,jend, * istart,iend) endif c Time series call puttimeseriesnc(ivartsid,infilenamets,ntseries,ncid, * maxts,ntime,nfirstyear,timeout,maxyear,iserrtsid,ntimeinterval, * icolumnts) instat=nf_close(ncid) if ( instat .ne. 0 ) then write(6,*) 'Problem closing file:',ncid,instat endif stop end SUBROUTINE GETINPUTINFNC(idim,jdim,xgrid,levels,dz, * kdimax,idz,nvars,maxvars,ivars,varname,istart,iend, * jstart,jend,ntime,itimestart,itimestart2,itimeinterval, * infilename, outfilename,ncvar,cunit,nunit,globalfile, * ftitle,ifsize,ifirstyear,itimeend,itimeend2, * varnamefull,ncvarf,ntseries,nfieldtypets, * varnamets,nspvts,cunitts,nsputs,filenamets, * varnamefullts,nspvfts,maxts,icolumnts,landmask,noland) parameter (idim0=360,jdim0=180,kdimax0=137) parameter (xdim=360.,xdim2=180.,ydim2=90.) parameter (icmax=10,icolumndefault=3) character*200 outfilename,infilename(maxvars),landmask character*80 varname,cunit,globalfile,varnamefull character*80 varnamets,cunitts,varnamefullts character*100 ctext character*200 filenamets(maxts) character*500 ftitle dimension idz(kdimax),dz(kdimax),dz0(kdimax) dimension ivars(maxvars),nfieldtypets(maxts) dimension icolumnts(maxts) dimension icomma(icmax,2) open(9,file='./woatoncfile.d',status='old') c write(6,*) 'file open' c Get ocean variable and units call clearstring(varname,80) read(9,'(a)') varname c write(6,*) 'varname',varname call findlastchar(ncvar,varname,80) call clearstring(varnamefull,80) read(9,'(a)') varnamefull call findlastchar(ncvarf,varnamefull,80) c write(6,*) 'xx',ncvarf,varnamefull call clearstring(cunit,80) read(9,'(a)') cunit call findlastchar(nunit,cunit,80) c Get number of variables to write out read(9,*) nvars do 30 n=1,nvars c write(6,*)'n',n c Read in variable types c 0 - analyzed mean c 1 - statitical mean c 2 - number of observations c 3 - standard deviations c 4 - standard error of the mean c 5 - observed minus analyzed c 6 - field minus annual c 7 - grid point c 8 - standard error of objectively analyzed value c 9 - ocean heat content anomaly read(9,*) ivars(n) c Get full input file names call clearstring(infilename(n),200) read(9,'(a)') infilename(n) c write(6,*) 'infile',infilename(n) 30 continue c Get output file name call clearstring(outfilename,200) read(9,'(a)') outfilename c write(6,*) 'out',outfilename c Get title call clearstring(ftitle,500) read(9,'(a)') ftitle c write(6,*) ftitle call findlastchar(ifsize,ftitle,500) c time interval (months),time start (0 for mean climatologies) c time end (0 for mean climatologies) c otherwise itimestart,itimeend are years, c itimestart2, itimeend2 are in units of time interval c ifirstyear is year from which to start counting times read(9,*) itimeinterval, itimestart, itimestart2, * itimeend,itimeend2,ifirstyear c write(6,*) 'ifirstyear',ifirstyear if ( itimeinterval .le. 12 ) then imult=12./itimeinterval else imult=1 endif ntime=(itimeend-itimestart)*imult+ * (itimeend2-itimestart2)+1 c Insert time in output file name if necessary if ( itimestart .eq. 0 ) then if ( itimeinterval .eq. 12 ) then itimeinsert=0 elseif ( itimeinterval .eq. 3 ) then itimeinsert=(13+itimestart2)-1 elseif ( itimeinterval .eq. 1 ) then itimeinsert=itimestart2 endif call insertnum('X',itimeinsert,outfilename,200) else call insertyear(itimestart,'!',outfilename,200) itimeendx=itimeend if ( itimeinterval .ge. 12 ) * itimeendx=itimeend+(itimeinterval/12)-1 call insertyear(itimeendx,'#',outfilename,200) itimestart2x=itimestart2 itimeend2x=itimeend2 if ( itimeinterval .lt. 12 ) then itimestart2x=((itimestart2-1)*itimeinterval)+1 itimeend2x=itimestart2x+itimeinterval-1 endif call insertnumexpand('?',itimestart2x,outfilename,200) call insertnumexpand('^',itimeend2x,outfilename,200) endif c Depths: which depth levels call clearstring(ctext,100) read(9,'(a)') ctext call findlastchar(nspd,ctext,100) call getwoadepthset(idz,dz,dz0,levels,ctext, * nspd,ilevitus,kdimax) c Starting/ending lat/lons, grid size read(9,*) xstart,xend,ystart,yend,xgrid c write(6,*) 'xstart',xstart xadd=xgrid/2. ystart=ystart+xadd if ( ystart .gt. ydim2 ) ystart=ydim2 yend=yend+xadd if ( yend .gt. ydim2 ) yend=ydim2 xstart=xstart+xadd if ( xstart .gt. xdim2 ) xstart=xstart-xdim xend=xend+xadd if ( xend .gt. xdim2 ) xend=xend-xdim call grid(ystart,xstart,jstart,istart,xgrid) call grid(yend,xend,jend,iend,xgrid) if ( istart .eq. iend ) iend=istart-1 idim=idim0/xgrid jdim=jdim0/xgrid c write(6,*) istart,iend,jstart,jend,xgrid c Get global attributes input file name call clearstring(globalfile,80) read(9,'(a)') globalfile c write(6,*) 'glob',globalfile c Get land mask file name call clearstring(landmask,200) read(9,'(a)') landmask read(9,*) noland c write(6,*) 'l',noland c Get time series info (if present) ntseries=0 read(9,*,end=4) ntseries call clearstring(varnamets,80) read(9,'(a)') varnamets call findlastchar(nspvts,varnamets,80) call clearstring(varnamefullts,80) read(9,'(a)') varnamefullts call findlastchar(nspvfts,varnamefullts,80) call clearstring(cunitts,80) read(9,'(a)') cunitts call findlastchar(nsputs,cunitts,80) do 60 nn=1,ntseries call clearstring(ctext,100) read(9,'(a)') ctext c write(6,*) ctext call findlastchar(nsp,ctext,100) call spacefindl(nsp,icmax,ncomma,icomma,ctext) c write(6,*) 'ncomma',ncomma read(ctext(icomma(1,1):icomma(1,1)+icomma(1,2)-1),*) * nfieldtypets(nn) c write(6,*) nfieldtypets(nn),icomma(2,1),icomma(2,2) if ( ncomma .gt. 1 ) then c write(6,*) 'ctext',ctext(icomma(2,1):icomma(2,1)+icomma(2,2)-1) read(ctext(icomma(2,1):icomma(2,1)+icomma(2,2)-1),*) * icolumnts(nn) else icolumnts(nn)=icolumndefault endif c write(6,*) 'here' call clearstring(filenamets(nn),200) read(9,'(a)') filenamets(nn) 60 continue 4 continue close(9) return end SUBROUTINE FILEOPENNC_CLASSIC(ncid,filename) C FILEOPENNC OPENS A NETCDF FILE FOR WRITING character*(*) filename character*120 systemname include 'netcdf.inc' write(6,*) 'filename:',filename c instat= nf_create(filename,NF_WRITE,ncid) c instat=nf_create(filename,NF_NETCDF4,ncid) instat=nf_create(filename,OR(NF_NETCDF4,NF_CLASSIC_MODEL),ncid) c overwrite if necessary if ( instat .eq. -35 ) then call clearstring(systemname,100) systemname='rm '//filename call system(systemname) c instat= nf_create(filename,NF_WRITE,ncid) instat=nf_create(filename,OR(NF_NETCDF4,NF_CLASSIC_MODEL),ncid) c instat=nf_create(filename,NF_NETCDF4,ncid) endif if ( instat .ne. 0 ) then write(6,*) 'Problem opening file:',filename write(6,*) instat stop endif return end SUBROUTINE SETUPFIELDVARNC(nfieldtype,varname,nspv,cunit, * nunit,ivarid,ivardims,ivardimsize,ncid,intfill,rfill, * intfile,varnamefull,nspvf,iallintegral,icheck) C SETUPFIELDVARNC SETS UP EACH TYPE OF VARIABLE FIELD C FOR NETCDF FILE parameter (maxfield=100) character*5 cfieldsuf(maxfield) character*15 firsttime character*50 cfieldcell(maxfield),cfieldname(maxfield) character*30 cfieldunit(maxfield) character*80 varname,outname,coutunit,cunit,varnamefull character*160 vfout character*500 ctext,cfieldcomm(maxfield) character*600 cout dimension ifield(maxfield),ncsuf(maxfield) dimension ncname(maxfield),nccomm(maxfield) dimension nccell(maxfield),ncunit(maxfield) dimension ifill(maxfield),ivardims(4),integral(maxfield) dimension icomma(maxfield,2),nsplit(maxfield) dimension ivardimsize(4) save firsttime,cfieldsuf,cfieldunit,cfieldcell,cfieldname save ifield,ncsuf,ncname,nccomm,nccell,ncunit,nfields,ifill, * integral,nsplit include 'netcdf.inc' if ( firsttime .ne. 'setupfieldvarnc' ) then firsttime='setupfieldvarnc' open(9,file='./woa_field_definitionsnc.dat',status='old') nfields=0 do 30 n=1,maxfield call clearstring(ctext,500) read(9,'(a)',end=44) ctext call findlastchar(nsp,ctext,500) call commafindx(nsp,maxfield,ncomma,icomma,ctext) c write(6,*) 'ncx',icomma(1,2),icomma(2,2),icomma(3,2) read(ctext(icomma(1,1):icomma(1,1)+icomma(1,2)-1),*) ifield(n) cfieldsuf(n)=ctext(icomma(2,1):icomma(2,1)+icomma(2,2)-1) ncsuf(n)=icomma(2,2) cfieldname(n)=ctext(icomma(3,1):icomma(3,1)+icomma(3,2)-1) ncname(n)=icomma(3,2) c write(6,*) 'ncname',n,ncname(n),ncomma,icomma(3,1),icomma(3,2) cfieldcomm(n)=ctext(icomma(4,1):icomma(4,1)+icomma(4,2)-1) nccomm(n)=icomma(4,2) nsplit(n)=-1 do 35 nn=1,nccomm(n)-3 35 if ( cfieldcomm(n)(nn:nn+3) .eq. 'VARX' ) nsplit(n)=nn c write(6,*) 'VARX',nsplit(n) if ( icomma(5,2) .gt. 0 .and. ncomma .ge. 5 ) then cfieldcell(n)=ctext(icomma(5,1):icomma(5,1)+icomma(5,2)-1) nccell(n)=icomma(5,2) else nccell(n)=0 endif if ( icomma(6,2) .gt. 0 .and. ncomma .ge. 6 ) then cfieldunit(n)=ctext(icomma(6,1):icomma(6,1)+icomma(6,2)-1) ncunit(n)=icomma(6,2) else ncunit(n)=0 endif if ( icomma(7,2) .gt. 0 .and. ncomma .ge. 7 ) then read(ctext(icomma(7,1):icomma(7,1)+icomma(7,2)-1),*) ifill(n) else ifill(n)=0 endif if ( icomma(8,2) .gt. 0 .and. ncomma .ge. 8 ) then read(ctext(icomma(8,1):icomma(8,1)+icomma(8,2)-1),*) * integral(n) else integral(n)=0 endif c write(6,*) 'xxxx',ncomma,icomma(8,2),integral(n) nfields=nfields+1 30 continue 44 continue close(9) endif intfile=0 ifieldfound=0 do 50 n=1,nfields if ( nfieldtype .eq. ifield(n) ) then if ( icheck .eq. 1 ) then if ( integral(n) .eq. 1 ) iallintegral=iallintegral+1 return endif c Variable definition call clearstring(outname,80) outname=varname(1:nspv)//"_"//cfieldsuf(n)(1:ncsuf(n)) nchar=nspv+ncsuf(n)+1 c Make sure variable name has no '-', as MATLAB dosnt like this do 77 nnx=1,nchar if ( outname(nnx:nnx) .eq. '-' ) outname(nnx:nnx)='_' 77 continue c write(6,*) 'zzz',outname(1:nchar) if ( ifill(n) .eq. 0 ) then instat=nf_def_var(ncid,outname(1:nchar),NF_FLOAT,4,ivardims, * ivarid) else instat=nf_def_var(ncid,outname(1:nchar),NF_INT,4,ivardims, * ivarid) endif if ( instat .ne. 0 ) then write(6,*) 'field variable problem',n, * instat endif c Set up chunking/compression - chunking will be entire variable field c (hence ivardimsize=chunksize). This is for later compression c instat=nf_def_var_chunking(ncid,ivarid,NF_CHUNKED,ivardimsize) c if ( instat .ne. 0 ) then c write(6,*) 'chunking definition problem',n, c * instat c endif c instat=nf_def_var_deflate(ncid,ivarid,0,1,2) c if ( instat .ne. 0 ) then c write(6,*) 'deflation setup problem',n, c * instat c endif c Standard name if ( integral(n) .eq. 0 .and. ifield(n) .ne. 3 ) then vfout=varnamefull(1:nspvf) nvout=nspvf if ( ifield(n) .eq. 2 .or. ifield(n) .eq. 4 ) then vfout=varnamefull(1:nspvf)//" "// * cfieldname(n)(1:ncname(n)) nvout=nspvf+ncname(n)+1 endif instat=nf_put_att_text(ncid,ivarid,'standard_name', * nvout,vfout(1:nvout)) if ( instat .ne. 0 ) then write(6,*) 'attribute standard name field problem', * instat endif endif c Long name instat=nf_put_att_text(ncid,ivarid,'long_name', * ncname(n),cfieldname(n)(1:ncname(n))) if ( instat .ne. 0 ) then write(6,*) 'attribute long name field problem', * instat endif c coordinates instat=nf_put_att_text(ncid,ivarid,'coordinates', * 18,'time lat lon depth' ) if ( instat .ne. 0 ) then write(6,*) 'attribute coordinate field problem', * instat endif c Comments cout=cfieldcomm(n)(1:nccomm(n)) ncout=nccomm(n) if ( nsplit(n) .gt. 0 ) then cout=cfieldcomm(n)(1:nsplit(n)-1)// * varnamefull(1:nspvf)// * cfieldcomm(n)(nsplit(n)+4:nccomm(n)) ncout=nccomm(n)-4+nspvf endif c instat=nf_put_att_text(ncid,ivarid,'comment', c * nccomm(n),cfieldcomm(n)(1:nccomm(n))) c write(6,*) 'cout',ncout,nccomm(n),nspvf,cout(1:ncout) c write(6,*) 'xxc',nsplit(n),cfieldcomm(n)(nsplit(n)+4:nccomm(n)) instat=nf_put_att_text(ncid,ivarid,'long_name', * ncout,cout) if ( instat .ne. 0 ) then write(6,*) 'attribute comment field problem', * instat endif c Cell method instat=nf_put_att_text(ncid,ivarid,'cell_methods', * nccell(n),cfieldcell(n)(1:nccell(n))) if ( instat .ne. 0 ) then write(6,*) 'attribute cell method field problem', * instat endif c Units if ( ncunit(n) .eq. 7 .and. cfieldunit(n)(1:7) .eq. * 'DEFAULT' ) then coutunit(1:nunit)=cunit(1:nunit) nout=nunit else coutunit(1:ncunit(n))=cfieldunit(n)(1:ncunit(n)) nout=ncunit(n) endif instat=nf_put_att_text(ncid,ivarid,'grid_mapping', * 3,'crs') if ( instat .ne. 0 ) then write(6,*) 'attribute unit field problem', * instat endif instat=nf_put_att_text(ncid,ivarid,'units', * nout,coutunit(1:nout)) if ( instat .ne. 0 ) then write(6,*) 'attribute unit field problem', * instat endif if ( ifill(n) .eq. 1 ) then instat=nf_put_att_int(ncid,ivarid,'_FillValue', * NF_INT,1,intfill) intfile=1 else instat=nf_put_att_real(ncid,ivarid,'_FillValue', * NF_REAL,1,rfill) intfile=0 endif if ( instat .ne. 0 ) then write(6,*) 'attribute fill value field problem', * instat endif ifieldfound=1 endif 50 continue if ( ifieldfound .eq. 0 ) then write(6,*) 'Field type not found:',nfieldtype endif return end SUBROUTINE LATITUDESETUP(ncid,jdimid,ilatid, * nbounddimid,ilatboundid) C SETS UP LATITUDE VARIABLE character*80 ctext dimension index(2) include 'netcdf.inc' instat=nf_def_var(ncid,'lat',NF_FLOAT,1,jdimid,ilatid) if ( instat .ne. 0 ) then write(6,*) 'latitude array problem',instat endif instat=nf_put_att_text(ncid,ilatid,'standard_name', * 8,'latitude') if ( instat .ne. 0 ) then write(6,*) 'attribute standard name latitude problem', * instat endif instat=nf_put_att_text(ncid,ilatid,'long_name', * 8,'latitude') if ( instat .ne. 0 ) then write(6,*) 'attribute long name latitude problem', * instat endif instat=nf_put_att_text(ncid,ilatid,'units', * 13,'degrees_north') if ( instat .ne. 0 ) then write(6,*) 'attribute units latitude problem', * instat endif instat=nf_put_att_text(ncid,ilatid,'axis', * 1,'Y') if ( instat .ne. 0 ) then write(6,*) 'attribute axis latitude problem', * instat endif instat=nf_put_att_text(ncid,ilatid,'bounds', * 8,'lat_bnds') if ( instat .ne. 0 ) then write(6,*) 'attribute bounds latitude problem', * instat endif index(1)=nbounddimid index(2)=jdimid instat=nf_def_var(ncid,'lat_bnds',NF_FLOAT,2, * index,ilatboundid) if ( instat .ne. 0 ) then write(6,*) 'latitudebounds array problem',instat endif call clearstring(ctext,80) ctext="latitude bounds" call findlastchar(nsp,ctext,80) instat=nf_put_att_text(ncid,ilatboundid,'comment', * nsp,ctext) if ( instat .ne. 0 ) then write(6,*) 'attribute comment lat bounds problem', * instat endif return end SUBROUTINE LONGITUDESETUP(ncid,idimid,ilonid,nbounddimid, * ilonboundid) C SETS UP LONGITUDE VARIABLE character*80 ctext dimension index(2) include 'netcdf.inc' instat=nf_def_var(ncid,'lon',NF_FLOAT,1,idimid,ilonid) if ( instat .ne. 0 ) then write(6,*) 'longitude array problem',instat endif instat=nf_put_att_text(ncid,ilonid,'standard_name', * 9,'longitude') if ( instat .ne. 0 ) then write(6,*) 'attribute standard name longitude problem', * instat endif instat=nf_put_att_text(ncid,ilonid,'long_name', * 9,'longitude') if ( instat .ne. 0 ) then write(6,*) 'attribute long name longitude problem', * instat endif instat=nf_put_att_text(ncid,ilonid,'units', * 12,'degrees_east') if ( instat .ne. 0 ) then write(6,*) 'attribute units longitude problem', * instat endif instat=nf_put_att_text(ncid,ilonid,'axis', * 1,'X') if ( instat .ne. 0 ) then write(6,*) 'attribute axis longitude problem', * instat endif instat=nf_put_att_text(ncid,ilonid,'bounds', * 8,'lon_bnds') if ( instat .ne. 0 ) then write(6,*) 'attribute bounds longitude problem', * instat endif index(1)=nbounddimid index(2)=idimid instat=nf_def_var(ncid,'lon_bnds',NF_FLOAT,2, * index,ilonboundid) if ( instat .ne. 0 ) then write(6,*) 'longitudebounds array problem',instat endif call clearstring(ctext,80) ctext="longitude bounds" call findlastchar(nsp,ctext,80) instat=nf_put_att_text(ncid,ilonboundid,'comment', * nsp,ctext) if ( instat .ne. 0 ) then write(6,*) 'attribute comment lon bounds problem', * instat endif return end SUBROUTINE TIMESETUP(ncid,itimedimid,itimeid, * nbounddimid,itimeboundid,ifirstyear) C SETS UP TIME VARIABLE character*32 cout character*80 ctext dimension index(2) include 'netcdf.inc' instat=nf_def_var(ncid,'time',NF_FLOAT,1,itimedimid,itimeid) if ( instat .ne. 0 ) then write(6,*) 'time array problem',instat endif instat=nf_put_att_text(ncid,itimeid,'standard_name', * 4,'time') if ( instat .ne. 0 ) then write(6,*) 'attribute standard name time problem', * instat endif instat=nf_put_att_text(ncid,itimeid,'long_name', * 4,'time') if ( instat .ne. 0 ) then write(6,*) 'attribute long name time problem', * instat endif cout='months since !!!!-01-01 00:00:00' call insertnum('!',ifirstyear,cout,32) instat=nf_put_att_text(ncid,itimeid,'units', * 32,cout) if ( instat .ne. 0 ) then write(6,*) 'attribute units time problem', * instat endif instat=nf_put_att_text(ncid,itimeid,'axis', * 1,'T') if ( instat .ne. 0 ) then write(6,*) 'attribute axis time problem', * instat endif instat=nf_put_att_text(ncid,itimeid,'climatology', * 18,'climatology_bounds') if ( instat .ne. 0 ) then write(6,*) 'attribute bounds time problem', * instat endif index(1)=nbounddimid index(2)=itimedimid instat=nf_def_var(ncid,'climatology_bounds',NF_FLOAT,2, * index,itimeboundid) if ( instat .ne. 0 ) then write(6,*) 'timebounds array problem',instat endif call clearstring(ctext,80) ctext="This variable defines the bounds of the"// * " climatological time period for each time step." call findlastchar(nsp,ctext,80) instat=nf_put_att_text(ncid,itimeboundid,'comment', * nsp,ctext) if ( instat .ne. 0 ) then write(6,*) 'attribute standard name time problem', * instat endif return end SUBROUTINE DEPTHSETUP(ncid,izdimid,izid, * nbounddimid,izboundid,integral) C SETS UP DEPTH VARIABLE character*80 ctext dimension index(2) include 'netcdf.inc' if ( integral .eq. 0 ) then instat=nf_def_var(ncid,'depth',NF_FLOAT,1,izdimid,izid) if ( instat .ne. 0 ) then write(6,*) 'depth array problem',instat endif instat=nf_put_att_text(ncid,izid,'standard_name', * 5,'depth') if ( instat .ne. 0 ) then write(6,*) 'attribute standard name depth problem', * instat endif instat=nf_put_att_text(ncid,izid,'bounds', * 10,'depth_bnds') if ( instat .ne. 0 ) then write(6,*) 'attribute bounds depth problem', * instat endif endif index(1)=nbounddimid index(2)=izdimid instat=nf_def_var(ncid,'depth_bnds',NF_FLOAT,2, * index,izboundid) if ( instat .ne. 0 ) then write(6,*) 'depthbounds array problem',instat endif call clearstring(ctext,80) ctext="depth bounds" call findlastchar(nsp,ctext,80) instat=nf_put_att_text(ncid,izboundid,'comment', * nsp,ctext) if ( instat .ne. 0 ) then write(6,*) 'attribute comment depth bounds problem', * instat endif if ( integral .eq. 1 ) izid=izboundid instat=nf_put_att_text(ncid,izid,'positive', * 4,'down') if ( instat .ne. 0 ) then write(6,*) 'attribute positive depth problem', * instat endif instat=nf_put_att_text(ncid,izid,'units', * 6,'meters') if ( instat .ne. 0 ) then write(6,*) 'attribute units depth problem', * instat endif instat=nf_put_att_text(ncid,izid,'axis', * 1,'Z') if ( instat .ne. 0 ) then write(6,*) 'attribute axis depth problem', * instat endif return end SUBROUTINE PUTDEPTHNC(ncid,levels,dz,kdimax,izid, * izboundid,integral,levels0) C PUTDEPTH ENTERS DEPTH INFORMATION INTO DEPTH AND C DEPTH BOUND ARRAY dimension index(2) dimension dz(kdimax) include 'netcdf.inc' do 50 k=1,levels0 if ( k .gt. 1 ) then dzb1=-(0.5*(dz(k)-dz(k-1)))+dz(k) else dzb1=dz(k) endif if ( k .lt. levels0 ) then dzb2=0.5*(dz(k+1)-dz(k))+dz(k) else dzb2=dz(levels) endif if ( integral .eq. 0 ) then instat=nf_put_var1_real(ncid,izid,k,dz(k)) if ( instat .ne. 0 ) then write(6,*) 'Problem entering depth:',dz(k),izid, * instat,ncid,k endif endif index(2)=k index(1)=1 instat=nf_put_var1_real(ncid,izboundid,index,dzb1) if ( instat .ne. 0 ) then write(6,*) 'Problem entering depth bound:',dz(k),izboundid, * instat,ncid,k,dzb1 endif index(1)=2 instat=nf_put_var1_real(ncid,izboundid,index,dzb2) if ( instat .ne. 0 ) then write(6,*) 'Problem entering depth bound:',dz(k),izboundid, * instat,ncid,k,dzb2 endif 50 continue return end SUBROUTINE PUTLONGITUDEnc(ncid,istart,iend,xgrid,ilonid, * ilonboundid) C PUTLONGITUDE ENTERS LONGITUDE INFORMATION INTO LONGITUDE AND C LONGITUDE BOUND ARRAY dimension index(2) include 'netcdf.inc' c write(6,*) 'lonbound3',ilonboundid j=180./xgrid idim=360./xgrid xgrid2=xgrid/2. iindex=0 ntime=1 if ( istart .gt. iend ) ntime=2 do 45 nt=1,ntime istart1=istart iend1=iend if ( ntime .eq. 2 ) then if ( nt .eq. 1 ) then istart1=istart iend1=idim else istart1=1 iend1=iend endif endif do 50 i=istart1,iend1 call ungrid(rlat,rlon,j,i,xgrid) rlonb1=rlon-xgrid2 rlonb2=rlon+xgrid2 iindex=iindex+1 instat=nf_put_var1_real(ncid,ilonid,iindex,rlon) if ( instat .ne. 0 ) then write(6,*) 'Problem entering longitude:',rlon,ilonid, * instat,ncid,iindex endif index(2)=iindex index(1)=1 instat= nf_put_var1_real(ncid,ilonboundid,index,rlonb1) if ( instat .ne. 0 ) then write(6,*) 'Problem entering longitude bound:',rlon, * ilonboundid,instat,ncid,iindex,rlonb1 endif index(1)=2 instat=nf_put_var1_real(ncid,ilonboundid,index,rlonb2) if ( instat .ne. 0 ) then write(6,*) 'Problem entering longitude bound:',rlon, * ilonboundid,instat,ncid,iindex,rlonb2 endif 50 continue 45 continue return end SUBROUTINE PUTLATITUDENC(ncid,jstart,jend,xgrid,ilatid, * ilatboundid) C PUTLATITUDE ENTERS LATITUDE INFORMATION INTO LATITUDE AND C LATITUDE BOUND ARRAY dimension index(2) include 'netcdf.inc' i=180./xgrid xgrid2=xgrid/2. do 50 j=jstart,jend call ungrid(rlat,rlon,j,i,xgrid) rlatb1=rlat-xgrid2 rlatb2=rlat+xgrid2 jindex=(j-jstart)+1 instat=nf_put_var1_real(ncid,ilatid,jindex,rlat) if ( instat .ne. 0 ) then write(6,*) 'Problem entering latitude:',rlat,ilatid, * instat,ncid,jindex endif index(2)=jindex index(1)=1 instat= nf_put_var1_real(ncid,ilatboundid,index,rlatb1) if ( instat .ne. 0 ) then write(6,*) 'Problem entering latitude bound:',rlat,ilatboundid, * instat,ncid,jindex,rlatb1 endif index(1)=2 instat=nf_put_var1_real(ncid,ilatboundid,index,rlatb2) if ( instat .ne. 0 ) then write(6,*) 'Problem entering latitude bound:',rlat,ilatboundid, * instat,ncid,jindex,rlatb2 endif 50 continue return end SUBROUTINE PUTTIMENC(ncid,ntime,ntimestart,ntimestart2, * ntimeinterval, itimeid,itimeboundid,nfirstyear,timeout, * maxyear) C PUTTIME ENTERS TIME INFORMATION INTO TIME AND C TIME BOUND ARRAY dimension index(2),timeout(maxyear) include 'netcdf.inc' xinterval2=ntimeinterval/2. xinterval2x=xinterval2 if ( ntimeinterval .ge. 12 ) * xinterval2=xinterval2/12. ntimeyear=12/ntimeinterval do 50 k=1,ntime ntimefull=0 ntimemod=0 if ( ntimeyear .gt. 1 ) then ntimefull=(ntimestart2+k-2)/ntimeyear ntimemod=((ntimestart2+k-2)- * (ntimefull*ntimeyear))*ntimeinterval c xtime=(ntimestart-nfirstyear)*12+(ntimefull*ntimeyear) xtime=(ntimestart-nfirstyear)*12+(ntimefull*12) * +ntimemod+xinterval2 else k1=ntimestart+k-1 xtime=(k1-nfirstyear)+xinterval2 c xtime=((k1-nfirstyear)*ntimeinterval)+xinterval2 c write(6,*) 'xtime',k1,ntimestart,k,nfirstyear, c * ntimeinterval,xinterval2 endif c write(6,*) 'xtime',ntimestart,ntimestart2,ntimeinterval, c * ntimefull,ntimemod,k,nfirstyear,ntimeyear,k c xtime=(ntimestart+k-1)+((ntimestart2-1)*ntimeinterval) c * +xinterval2 timeout(k)=xtime if ( ntimeyear .le. 1 ) then xtime=xtime*12. endif xtimeb1=xtime-xinterval2x xtimeb2=xtime+xinterval2x instat=nf_put_var1_real(ncid,itimeid,k,xtime) if ( instat .ne. 0 ) then write(6,*) 'Problem entering time:',xtime,itimeid, * instat,ncid,k endif index(2)=k index(1)=1 instat=nf_put_var1_real(ncid,itimeboundid,index,xtimeb1) if ( instat .ne. 0 ) then write(6,*) 'Problem entering time bound:',xtime,itimeboundid, * instat,ncid,k,xtimeb1 endif index(1)=2 instat=nf_put_var1_real(ncid,itimeboundid,index,xtimeb2) if ( instat .ne. 0 ) then write(6,*) 'Problem entering time bound:',xtime,itimeboundid, * instat,ncid,k,xtimeb2 endif 50 continue return end SUBROUTINE SETUP_GLOBALATT_WOANC(ncid,globalfile,jstart, * jend,istart,iend,xgrid,dz1,dz2,ftitle,ifsize, * itimestart,itimestart2,itimeinterval,itimeend, * itimeend2,outfilename) C SETUP_GLOBALATT_WOANC SETS UP STANDARD GLOBAL ATTRIBUTES C FOR WOA NETCDF FORM parameter (maxglobe=500,icmax=20) character*80 cname,globalfile character*500 ctext,cout character*(*) ftitle character*(*) outfilename dimension icomma(icmax,2) dimension idt(3) include 'netcdf.inc' call idate(idt) iday=idt(1) month=idt(2) iyear=idt(3) open(9,file=globalfile,status='old') iout=0 xout=0.0 jnsp=0 insp=0 xgrid2=xgrid/2. do 50 n=1,maxglobe call clearstring(ctext,500) write(6,*) 'going in',n read(9,'(a)',end=4) ctext write(6,*) ctext,n call findlastchar(nsp,ctext,500) call charfindl(nsp,icmax,ncomma,icomma,';',ctext) call clearstring(cname,80) cname=ctext(icomma(1,1):icomma(1,1)+icomma(1,2)-1) insp=icomma(1,2) read(ctext(icomma(2,1):icomma(2,1)+icomma(2,2)-1),*) itype call clearstring(cout,500) if ( ncomma .ge. 3 .and. icomma(3,2) .gt. 0 ) then cout=ctext(icomma(3,1):icomma(3,1)+icomma(3,2)-1) jnsp=icomma(3,2) if ( (jnsp .eq. 7 .and. cout(1:7) .eq. 'SPECIAL') .or. * (jnsp .gt. 7 .and. cout(1:7) .eq. 'SPECIAL' .and. * insp .eq. 2 .and. cname(1:insp) .eq. 'id') ) then if ( insp .eq. 5 .and. * cname(1:insp) .eq. 'title' ) then c open(10,file='nc_title.dat',status='old') call clearstring(cout,500) c read(10,'(a)') cout c call findlastchar(jnsp,cout,500) cout=ftitle(1:ifsize) jnsp=ifsize elseif ( insp .eq. 2 .and. * cname(1:insp) .eq. 'id' ) then call findlastchar(nspid,outfilename,200) call charfindl(nspid,icmax,ncomma,icomma,'/',outfilename) c write(6,*) nspid,'|',outfilename(1:nspid),'|',ncomma nspfo=jnsp-7 if ( nspfo .gt. 0 ) then c write(6,*) 'nspfo',nspfo,ncomma,icomma(ncomma,1), c * icomma(ncomma,2),jnsp,cout(1:jnsp) cout= * outfilename(icomma(ncomma,1): * icomma(ncomma,1)+icomma(ncomma,2)-1)// * cout(8:jnsp) jnsp=icomma(ncomma,2)+nspfo else cout= * outfilename(icomma(ncomma,1): * icomma(ncomma,1)+icomma(ncomma,2)-1) jnsp=icomma(ncomma,2) endif elseif ( insp .eq. 18 .and. * cname(1:insp) .eq. 'geospatial_lat_min' ) then call ungrid(xout,rlon,jstart,istart,xgrid) xout=xout-xgrid2 elseif ( insp .eq. 18 .and. * cname(1:insp) .eq. 'geospatial_lat_max' ) then call ungrid(xout,rlon,jend,istart,xgrid) xout=xout+xgrid2 elseif ( insp .eq. 18 .and. * cname(1:insp) .eq. 'geospatial_lon_min' ) then call ungrid(rlat,xout,jend,istart,xgrid) xout=xout-xgrid2 elseif ( insp .eq. 18 .and. * cname(1:insp) .eq. 'geospatial_lon_max' ) then call ungrid(rlat,xout,jend,iend,xgrid) xout=xout+xgrid2 elseif ( insp .eq. 25 .and. ( * cname(1:insp) .eq. 'geospatial_lat_resolution' * .or. cname(1:insp) .eq. * 'geospatial_lon_resolution') ) then call clearstring(cout,80) write(cout(1:),800) xgrid 800 format(f4.2," degrees") jnsp=12 elseif ( insp .eq. 23 .and. * cname(1:insp) .eq. 'geospatial_vertical_min' ) * then xout=dz1 elseif ( insp .eq. 23 .and. * cname(1:insp) .eq. 'geospatial_vertical_max' ) * then xout=dz2 elseif ( (insp .eq. 12 .and. * cname(1:insp) .eq. 'date_created') .or. * (insp .eq. 13 .and. * cname(1:insp) .eq. 'date_modified') ) then cout='XXXX-YY-ZZ' call insertnum('X',iyear,cout,11) call insertnum('Y',month,cout,11) call insertnum('Z',iday,cout,11) jnsp=11 elseif ( insp .eq. 22 .and. * cname(1:insp) .eq. 'time_coverage_duration' .and. * itimestart .ne. 0 ) then itime1d=0 if ( itimestart2 .eq. 0 ) then itime1d=(itimeend+(itimeinterval/12)-1)-itimestart cout='P!!Y' jnsp=4 call insertnum('!',itime1d+1,cout,4) else itime1d=((itimestart2-1)*itimeinterval)+1 itime2d=itimeend2*itimeinterval if ( itimeend .eq. itimestart ) then itimespan=itime2d-itime1d+1 else itimespan=(((itimeend)*12-(12-itime2d))- * ((itimestart-1)*12+itime1d))+1 endif if ( itimespan .gt. 999 ) then cout='P!!!!M' jnsp=6 elseif ( itimespan .gt. 99 ) then cout='P!!!M' jnsp=5 else cout='P!!M' jnsp=4 endif c write(6,*) 'itimespan1',itimespan,cout,jnsp call insertnum('!',itimespan,cout,jnsp) c itimeyear=12/itimeinterval c if ( itime1d .gt. 0 ) then c itime2d=(itimeend2+(itimeyear-itimestart2)) c itime1d=(itime1d-1)*12 c else c itime2d=itimeend2-itimestart2 c endif c itime2d=(itime2d+1)*itimeinterval c cout='P!!M' c write(6,*) 'here',itime1d,itime2d c call insertnum('!',itime1d+itime2d,cout,4) endif c jnsp=4 elseif ( (insp .eq. 22 .and. * cname(1:insp) .eq. 'time_coverage_duration' .and. * itimestart .eq. 0 ) * .or. (insp .eq. 24 .and. * cname(1:insp) .eq. 'time_coverage_resolution') * ) then inotfullyear=0 if ( itimeinterval .ge. 12 ) then cout='P!!Y' it2=itimeinterval/12 if ( it2*12 .eq. itimeinterval ) then call insertnum('!',it2,cout,4) jnsp=4 else inotfullyear=1 endif endif if ( itimeinterval .lt. 12 .or. * inotfullyear .eq. 1 ) then cout='P!!M' c write(6,*) 'itimespan2',itimeinterval,cout,jnsp call insertnum('!',itimeinterval,cout,4) jnsp=4 endif elseif ( insp .eq. 19 .and. * cname(1:insp) .eq. 'time_coverage_start') then cout='!!!!-##-01' jnsp=10 if ( itimeinterval .ge. 12 ) then nin2=1 else nin2=((itimestart2-1)*itimeinterval)+1 endif call insertnum('!',itimestart,cout,jnsp) call insertnum('#',nin2,cout,jnsp) endif elseif ( itype .eq. 0 ) then read(ctext(icomma(3,1):icomma(3,1)+icomma(3,2)-1),*) iout elseif ( itype .eq. -1 ) then read(ctext(icomma(3,1):icomma(3,1)+icomma(3,2)-1),*) xout endif else jnsp=0 endif if ( itype .eq. 1 ) then instat=nf_put_att_text(ncid,NF_GLOBAL,cname(1:insp), * jnsp,cout(1:jnsp)) c write(6,*) 'instat',instat,cname(1:insp),insp,jnsp,cout(1:jnsp) elseif ( itype .eq. 0 ) then instat= nf_put_att_int(ncid,NF_GLOBAL,cname(1:insp), * NF_INT,1,iout) elseif ( itype .eq. -1 ) then instat= nf_put_att_real(ncid,NF_GLOBAL,cname(1:insp), * NF_FLOAT,1,xout) endif if ( instat .ne. 0 ) then write(6,*) 'Problem assigning global att',instat, * cname(1:insp),xout endif 50 continue 4 continue c write(6,*) 'done' close(9) return end SUBROUTINE SETWOAFILENAME(filename,itimestart,itimeinterval, * itimestart2,filename0,itz) C SETWOAFILENAME SETS WOA FILE NAME BASED ON TIME PERIOD C INPUT character*100 filename0,filename call clearstring(filename,100) filename=filename0 call findlastchar(nsp,filename,100) if ( itimestart .eq. 0 ) then if ( itimeinterval .eq. 12 ) then itimeinsert=0 elseif ( itimeinterval .eq. 3 ) then itimeinsert=(13+itimestart2+itz)-2 elseif ( itimeinterval .eq. 1 ) then itimeinsert=(itimestart2+itz)-1 endif call insertnum('X',itimeinsert,filename,nsp) c write(6,*) 'insert',itimeinsert,itimestart,itimeinterval else iyearadd=itz-1 if ( itimeinterval .lt. 12 ) then inyear=12/itimeinterval iyearadd=((itz+itimestart2)-2)/inyear im1=((((itz+itimestart2)-2)-(iyearadd*inyear))*itimeinterval)+1 im2=im1+(itimeinterval-1) c write(6,*) 'insert',im1,im2,itimeinterval,itimestart2,itz, c * iyearadd,inyear call insertnumexpand('?',im1,filename,100) call insertnumexpand('^',im2,filename,100) endif iyear1=itimestart+iyearadd iyear2=iyear1 if ( itimeinterval .ge. 12 ) then iyear2=iyear1+((itimeinterval/12)-1) endif call insertyear(iyear1,'!',filename,nsp) call insertyear(iyear2,'#',filename,nsp) if ( itimeinterval .lt. 12 ) then month1=mod((itz+itimestart2)-1,12/itimeinterval)+1 month2=month1+itimeinterval-1 call insertnumexpand('?',month1,filename,nsp) call insertnumexpand('^',month2,filename,nsp) endif endif call addCend(filename,100,1) return end SUBROUTINE READTIMESERIES(filename,xyear,xval,xse, * maxyear,nyear,icolumnts) C READTIMESERIES READS IN TIME SERIES DATA, INCLUDING C DATE (DECIMAL YEAR), INTEGRAL/MEAN VALUE, STANDARD C ERROR OF THE MEAN parameter (icmax=20) character*(*) filename character*100 ctext dimension xyear(maxyear),xval(maxyear),xse(maxyear) dimension icomma(icmax,2) open(10,file=filename,status='old') xmin=-3000. xmax=3000. nyear=0 do 50 nn=1,maxyear call clearstring(ctext,100) read(10,'(a)',end=4) ctext c write(6,*) ctext,icolumnts call findlastchar(nsp,ctext,100) call spacefindl(nsp,icmax,ncomma,icomma,ctext) read(ctext(icomma(1,1):icomma(1,1)+icomma(1,2)-1),*,err=50) * xyear(nn) c write(6,*) 'here',xyear(nn) ii=icolumnts read(ctext(icomma(ii,1):icomma(ii,1)+icomma(ii,2)-1),*) * xval(nn) ii=ii+1 read(ctext(icomma(ii,1):icomma(ii,1)+icomma(ii,2)-1),*) * xse(nn) c read(10,*,end=4) xyear(nn),xdum,xval(nn),xse(nn) nyear=nyear+1 50 continue 4 continue close(10) return end SUBROUTINE PUTTIMESERIESNC(ivartsid,infilenamets,ntseries,ncid, * maxts,ntime,nfirstyear,timeout,maxyear,iserrtsid,itimeinterval, * icolumnts) C PUTTIMESERIESNC PUTS TIME SERIES NETCDF FORM parameter (maxtime=12000) character*(*) infilenamets(maxts) dimension timeout(maxyear),xyear(maxtime) dimension xval(maxtime),xse(maxtime) dimension ivartsid(maxts),iserrtsid(maxts) dimension icolumnts(maxts) include 'netcdf.inc' do 30 nf=1,ntseries call readtimeseries(infilenamets(nf),xyear,xval,xse, * maxtime,nyear,icolumnts(nf)) do 50 n=1,nyear xy=xyear(n)-nfirstyear if ( itimeinterval .lt. 12 ) then xy=xy*12. c iyear=xyear(n) c xy0=xyear(n)-iyear c xy0=xy0*12. c xy=iyear+xy0 endif c write(6,*) 'fff',nyear,xyear(n),xy,xval(n),xse(n), c * ntime,timeout(1),timeout(ntime),nfirstyear,itimeinterval if ( xy .lt. timeout(1) .or. xy .gt. timeout(ntime) ) * goto 50 do 55 n2=1,ntime if ( xyear(n) .gt. 2012. ) then c write(6,*) 'xy',n2,timeout(n2) endif if ( xy .ge. timeout(n2) - 0.001 .and. * xy .le. timeout(n2)+0.001 ) then instat=nf_put_var1_real(ncid,ivartsid(nf),n2,xval(n)) if ( instat .ne. 0 ) then write(6,*) 'Problem entering time series:', * xval(n),ivartsid(nf), instat,ncid,nf,n endif instat=nf_put_var1_real(ncid,iserrtsid(nf),n2,xse(n)) if ( instat .ne. 0 ) then write(6,*) 'Problem entering time series se:', * xval(n),iserrtsid(nf), instat,ncid,nf,n endif endif 55 continue 50 continue 30 continue return end SUBROUTINE SETUPTIMESERIESNC(ntseries,nfieldtype, * cfieldnamets,ncfield,cunitts,ncunit,varname,nspv, * ivartsid,iserrtsid,maxts,itimedimid,ncid,idimid, * jdimid,ibasinid) C SETUPFIELDVARNC SETS UP EACH TYPE OF VARIABLE FIELD C FOR NETCDF FILE parameter (maxfield=100) character*20 firsttime character*2 cabbrev(maxfield) character*80 careaname(maxfield) character*500 ctext character*80 varname,outname character*80 cfieldnamets character*80 cunitts dimension ifield(maxfield) dimension icomma(maxfield,2) dimension ncname(maxfield),nfieldtype(maxts) dimension ivartsid(maxts),iserrtsid(maxts) dimension index(2) save firsttime,cabbrev,careaname save nfields,ifield,ncname include 'netcdf.inc' if ( firsttime .ne. 'setuptimeseriesvarnc' ) then firsttime='setuptimeseriesvarnc' open(9,file='./timeseries_area_info.txt', * status='old') nfields=0 do 30 n=1,maxfield call clearstring(ctext,500) read(9,'(a)',end=44) ctext call findlastchar(nsp,ctext,500) call commafindx(nsp,maxfield,ncomma,icomma,ctext) read(ctext(icomma(1,1):icomma(1,1)+icomma(1,2)-1),*) ifield(n) careaname(n)=ctext(icomma(2,1):icomma(2,1)+icomma(2,2)-1) ncname(n)=icomma(2,2) cabbrev(n)=ctext(icomma(3,1):icomma(3,1)+icomma(3,2)-1) nfields=nfields+1 30 continue 44 continue close(9) endif do 50 n=1,ntseries ifieldfound=0 do 55 n2=1,nfields if ( nfieldtype(n) .eq. ifield(n2) ) then c Variable definition call clearstring(outname,80) outname=varname(1:nspv)//"_"//cabbrev(n2) nchar=nspv+3 do 155 nn2=1,nchar 155 if ( outname(nn2:nn2) .eq. '-' ) outname(nn2:nn2)='_' if ( outname(1:4) .eq. '3mon' ) outname(1:4)='seas' if ( outname(1:5) .eq. '1yrav' ) outname(1:5)='yearl' write(6,*) 'zz',ncid,itimedimid,outname(1:nchar) instat=nf_def_var(ncid,outname(1:nchar),NF_FLOAT,1,itimedimid, * ivartsid(n)) if ( instat .ne. 0 ) then write(6,*) 'field variable problem',n, * instat endif call clearstring(outname,80) outname=varname(1:nspv)//"_se_"//cabbrev(n2) nchar=nspv+6 do 156 nn2=1,nchar 156 if ( outname(nn2:nn2) .eq. '-' ) outname(nn2:nn2)='_' if ( outname(1:4) .eq. '3mon' ) outname(1:4)='seas' if ( outname(1:5) .eq. '1yrav' ) outname(1:5)='yearl' instat=nf_def_var(ncid,outname(1:nchar),NF_FLOAT,1,itimedimid, * iserrtsid(n)) if ( instat .ne. 0 ) then write(6,*) 'field variable problem se',n, * instat endif c Standard name : no standard name c instat=nf_put_att_text(ncid,ivarid,'standard_name', c * nspvf,varnamefull(1:nspvf)) c if ( instat .ne. 0 ) then c write(6,*) 'attribute long name field problem', c * instat c endif c endif c Long name call clearstring(outname,80) c outname=careaname(n2)(1:ncname(n2))//"_"// outname=careaname(n2)(1:ncname(n2))// * cfieldnamets(1:ncfield) nchar=ncname(n)+ncfield+1 instat=nf_put_att_text(ncid,ivartsid(n),'long_name', * nchar,outname(1:nchar)) if ( instat .ne. 0 ) then write(6,*) 'attribute long name field problem', * instat endif call clearstring(outname,80) outname=careaname(n2)(1:ncname(n2))//"_standarderror"// * cfieldnamets(1:ncfield) nchar=ncname(n2)+ncfield+14 instat=nf_put_att_text(ncid,iserrtsid(n),'long_name', * nchar,outname(1:nchar)) if ( instat .ne. 0 ) then write(6,*) 'attribute long name field problem', * instat endif c Units instat=nf_put_att_text(ncid,ivartsid(n),'units', * ncunit,cunitts(1:ncunit)) if ( instat .ne. 0 ) then write(6,*) 'attribute unit field problem', * instat endif instat=nf_put_att_text(ncid,iserrtsid(n),'units', * ncunit,cunitts(1:ncunit)) if ( instat .ne. 0 ) then write(6,*) 'attribute unit field problem', * instat endif ifieldfound=1 endif 55 continue if ( ifieldfound .eq. 0 ) then write(6,*) 'time series type not found:',nfieldtype endif 50 continue c Basin mask if ( ntseries .gt. 0 ) then index(1)=idimid index(2)=jdimid instat=nf_def_var(ncid,'basin_mask',NF_INT,2,index, * ibasinid) if ( instat .ne. 0 ) then write(6,*) 'basin mask problem',n, * instat endif c Long name instat=nf_put_att_text(ncid,ibasinid,'long_name', * 10,'basin mask') if ( instat .ne. 0 ) then write(6,*) 'basin attribute long name field problem', * instat endif c Fill value instat=nf_put_att_int(ncid,ibasinid,'_FillValue', * NF_INT,1,-100) if ( instat .ne. 0 ) then write(6,*) 'basin attribute fill value field problem', * instat endif c Grid mapping instat=nf_put_att_text(ncid,ibasinid,'grid_mapping', * 3,'crs') if ( instat .ne. 0 ) then write(6,*) 'basin attribute grid mapping field problem', * instat endif c flag values call clearstring(outname,80) outname='-1, 1, 2, 3' call findlastchar(nsp,outname,80) instat=nf_put_att_text(ncid,ibasinid,'flag_values', * nsp,outname) if ( instat .ne. 0 ) then write(6,*) 'basin attribute flag values field problem', * instat endif c flag meaning call clearstring(outname,80) outname="Not_used_in_Basin_Integrals "// * "Atlantic_Basin Pacific_Basin Indian_Basin" call findlastchar(nsp,outname,80) instat=nf_put_att_text(ncid,ibasinid,'flag_values', * nsp,outname) if ( instat .ne. 0 ) then write(6,*) 'basin attribute flag meaning field problem', * instat endif c comments call clearstring(outname,80) outname="For North Hem. use Lats 0.5 to 89.5"// * " For South use Lats -89.5 to -0.5" call findlastchar(nsp,outname,80) instat=nf_put_att_text(ncid,ibasinid,'comment', * nsp,outname) if ( instat .ne. 0 ) then write(6,*) 'basin attribute flag meaning field problem', * instat endif endif return end SUBROUTINE PUTBASINDEFS(ncid,ibasinid,jstart,jend, * istart,iend) C GENERATE AND ENTER INTO NETCDF FILE BASIN DEFINITIONS parameter (idim=360,jdim=180) character*80 filename dimension iland(idim,jdim),ibasin(idim,jdim) dimension ivarindex(2) iatl=0 ipac=0 iind=0 inotused=0 ionland=0 ilandmask=-1 call extraname('masks'//CHAR(0),'etop1.d'//CHAR(0),filename) call fileopen(filename,'rb+'//CHAR(0),ilandmask) ibasinmask=-1 call extraname('masks'//CHAR(0),'basinmasknosouth.d'//CHAR(0), * filename) call fileopen(filename,'rb+'//CHAR(0),ibasinmask) call readleveli(1,1,1,1,1,idim,jdim,idim,jdim,idim,jdim, * 1,1,ilandmask,iland) call onefileclose(ilandmask) call readleveli(1,1,1,1,1,idim,jdim,idim,jdim,idim,jdim, * 1,1,ibasinmask,ibasin) call onefileclose(ibasinmask) write(6,*) 'iend',jstart,jend,istart,iend do 60 j=jstart,jend ivarindex(2)=(j-jstart+1) itimex=1 if ( istart .gt. iend ) itimex=2 index1=0 do 62 itx=1,itimex istart1=istart iend1=iend if ( itimex .gt. 1 ) then if ( itx .eq. 1 ) iend1=idim if ( itx .eq. 2 ) istart1=1 endif do 65 i=istart1,iend1 index1=index1+1 ivarindex(1)=index1 if ( iland(i,j) .le. 1 ) then iout=-100 ionland=ionland+1 elseif ( ibasin(i,j) .eq. 1 .or. ibasin(i,j) .eq. 11 ) then iout=1 iatl=iatl+1 elseif ( ibasin(i,j) .eq. 2 .or. ibasin(i,j) .eq. 12 ) then iout=2 ipac=ipac+1 elseif ( ibasin(i,j) .eq. 3 .or. ibasin(i,j) .eq. 56 ) then iout=3 iind=iind+1 else iout=-1 inotused=inotused+1 endif if ( iout .gt. -100 ) then instat=nf_put_var1_int(ncid,ibasinid,ivarindex,iout) if ( instat .ne. 0 ) then write(6,*) 'Problem writing out value',instat,iout, * varindex,i,j endif endif 65 continue 62 continue 60 continue write(6,*) 'ocean',ionland,iatl,ipac,iind,inotused return end SUBROUTINE SETUPCRSNC(ncid,icrsid,iprofileid) C SETUPCRSNC SETS UP LAT/LON INFORMATION SET include 'netcdf.inc' instat= nf_def_var(ncid,'crs',NF_INT,0,0, * icrsid) if ( instat .ne. 0 ) then write(6,*) 'Problem setting up crs',instat stop endif instat= nf_put_att_text(ncid,icrsid,'grid_mapping_name', * 18,'latitude_longitude') instat= nf_put_att_text(ncid,icrsid,'epsg_code', * 9,'EPSG:4326') instat= nf_put_att_real(ncid,icrsid, * 'longitude_of_prime_meridian',NF_FLOAT,1,0.0) instat= nf_put_att_real(ncid,icrsid, * 'semi_major_axis',NF_FLOAT,1,6378137.0) instat= nf_put_att_real(ncid,icrsid, * 'inverse_flattening',NF_FLOAT,1,298.257223563) c instat= nf_def_var(ncid,'profile',NF_INT,0,0, c * iprofiledim) if ( instat .ne. 0 ) then write(6,*) 'Problem setting up profile',instat stop endif c instat= nf_put_att_text(ncid,iprofileid,'cf_role', c * 10,'profile_id') c Global attributes c instat= nf_put_att_text(ncid,0,'standard_name_vocabulary', c * 6,'CF-1.5') c instat= nf_put_att_text(ncid,0,'featureType', c * 7,'profile') c instat= nf_put_att_text(ncid,0,'cdm_data_type', c * 7,'Profile') c instat= nf_put_att_text(ncid,0,'Conventions', c * 6,'CF-1.5') c instat= nf_put_att_text(ncid,0,'Metadata_Conventions', c * 30,'Unidata Dataset Discovery v1.0') return end SUBROUTINE GETWOADEPTHSET(idz,dz,dz0,levels,ctext, * nspd,ilevitus,kdimax) c GETWOADEPTHSET SETS DEPTHS AND INDICES TO BE C WRITTEN OUT TO FILE. C IF THERE IS AN 'L' IN THE FIRST SPACE, C GETS DEPTH SET INDICES FOR LEVITUS DEPTHS C WITHIN THE FULL DEPTH SET parameter (kdim=40,icmax=5) character*80 filename character*100 text character*(*) ctext dimension level2(kdim) dimension icomma(icmax,2) dimension idz(kdimax),dz(kdimax),dz0(kdimax) ilevitus=0 if ( ctext(1:1) .eq. 'L' ) then do 30 nn=1,nspd-1 30 ctext(nn:nn)=ctext(nn+1:nn+1) nspd=nspd-1 ilevitus=1 endif call spacefindl(nspd,icmax,ncomma,icomma,ctext) read(ctext(icomma(1,1):icomma(1,1)+icomma(1,2)-1),*) * indz1 read(ctext(icomma(2,1):icomma(2,1)+icomma(2,2)-1),*) * indz2 if ( ilevitus .eq. 0 ) then do 40 k=indz1,indz2 level2(k)=k 40 continue else c subset to levitus depths call clearstring(filename,80) call extraname("sys.inf"//CHAR(0), * "standepths_levtofull.txt"//CHAR(0),filename) open(11,file=filename,status='old') read(11,'(a)') text read(11,'(a)') text indexmax=0 do 45 kk=1,kdim call clearstring(text,80) read(11,'(a)',end=46) text call findlastchar(nspt,text,80) call commafindx(nspt,icmax,ncomma,icomma,text) read(text(icomma(2,1):icomma(2,1)+icomma(2,2)-1),*) * index1 read(text(icomma(3,1):icomma(3,1)+icomma(3,2)-1),*) * level2(index1) indexmax=index1 45 continue 46 continue if ( indexmax .lt. indz2 ) indz2=indexmax close(11) endif levels=(indz2-indz1)+1 call setstandepsfromfile(nztot,dz0) c write(6,*) 'depthin',nztot,dz0(levels) do 52 k=indz1,indz2 idz(k-indz1+1)=level2(k) dz(k-indz1+1)=dz0(level2(k)) 52 continue c write(6,*) 'dz',dz(levels) return end