program basin_layer10ann_climatology c@***h* OCEAN_HEAT_CONTENT/basin_layer10ann_climatology.f c c NAME c basin_layer10ann_climatology.f c c PURPOSE c Specifically for OCEAN_HEAT_CONTENT (OHC), OHC calculated c for each grid box over designated ocean layer (in program c layer_field_10ann_climatology.f) and standard error of OHC c over same (in program layer_fieldSE_K05.f) are summed over c designated latitude/longitude range globally as well as for c the Atlantic, Pacific, and Indian Ocean basins to give an c integrated OHC. c c DESCRIPTION c basin_layer10ann_climatology.f is a multipurpose program which c either sums or averages over given latitude/longitude range an c ocean variable read in from input file. c c INPUT c The main input to basin_layer10ann_climatology.f is the OHC over c designated layer calculated in layer_field_10ann_climatology.f and c the standard deviation calculated in layer_fieldSE_K05.f. User input, c such as form of input file, whether to itnegrate over average, time c span (monthly, quarterly, yearly, pentadal), time period c (specific time span), latitude/longitude range, name c of output file are found in basin_layer10ann-climatology.inf. c Designated basins to use (1) or not use (0) c are listed in OCL_basins.inf. at [MAINDIR]/sys.inf/areamask.1deg c contains definitions (lat/lon bounds) for each defined basin. c MAINDIR is the main directory for the OHC calculations under c which all subdirectories are found. c c OUTPUT c For OHC, OHC integrated over designated latitude/longitude c region for given basins for given time span, over designated c time period as well as standard error of the mean, in an c ASCII file. c cAUTHOR c John Antonov for 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 8, 2017 Initial version prepared for Reference Environmental c Data Record (REDR) program c c***** C 02.28.07 updated to work with monthly anomalies C 01/20/06 extended list of OCL basins (up to 100) C 01/09/06 using .se and .2se input files C 10/26/04: for AO, IO, PO, WO or User Basin C for each LAYER and ocean BASIN C one output file with timeseries of C WEIGHTED/UNWEIGHTED, TRIMMED/UNTRIMMED C BASIN TOTAL or MEAN VALUES C C 08/11/04 updated SE for basin total C 03/29/04 for monthly climatology output time:1. not 1901.5 C new parameter KID: C record number for global field in file C C Basins are defined according to OCL basins mask and/or C latitudinal and/or longitudinal bounds C C 06/12/01 C C TorM=1 basin total C = basin mean C C WorUW=1 weighted basin value C =2 unweighted basin value C C AorV=0 weights equal 1. (volume files would not be opened) C =1 weights are ODSQ areas C =2 weights are ODSQ volumes of layer C C SEtype=1 (se of basin value is sum(w*ODSQse) C 2 (se of basin value is sqrt(sum((w*ODSQse)**2))/sum(w) C 0 none C C 04/15/00 C C a0 area of 1x1 at equator in 10**10 m**2 C parameter (kid=1,kdim2=137) parameter (idim=360, jdim=180, kdim=102) parameter (bmiss=-1.E10, nrec=idim*jdim) parameter (pi=3.14159265, pirat=pi/180.) parameter (hmsph=360.0, h05=0.5, h1=1.) parameter (a0=1.237) parameter (icmax=5) dimension icomma(icmax,2) dimension f(idim,jdim), f1(idim,jdim), f2(idim,jdim), arw(jdim,4) dimension itopo(idim,jdim) dimension trim(nrec), wtrim(nrec), y(nrec), yw(nrec), ferr(nrec) dimension fvol(idim,jdim), volk(idim,jdim) dimension fmask(idim,jdim,kdim) dimension iup(30),ideep(30),glub(kdim),hall(kdim2) character*200 infl, inse, ingrand, out, mess, involyr, involkl character*200 junk, bOCL, umask(5) character*4 cdum, namend(3) character*3 ima, csouth, cnorth character*2 ctemp, cyr, cyr5 character*1 ox(5) character*6 dotend,doten1 character*4 cam(kdim) character*6 cstr(10) character*5 c6mon(23) character*200 ctext integer torm,woruw,aorv,setype trik(au)=au+hmsph*int(h05*(h1-sign(h1,au))) odsqa(alat)=abs(cos(alat)) data f/nrec*bmiss/, f1/nrec*bmiss/, f2/nrec*bmiss/ data ferr/nrec*0./ data zh/0./,azh/0./,nz/0/, se/0./, vol/0./ data fvol/nrec*0./,volk/nrec*0./, areabzn/0./ data c6mon/' ', * '1-6','7-12', * '1-3','4-6','7-9','10-12', * ' ',' ',' ',' ', * '1-1','2-2','3-3','4-4','5-5','6-6','7-7','8-8', * '9-9','10-10','11-11','12-12'/ data ox/'a','i','p','w','U'/ c Get depths call setstandepsfromfile(numdeps,hall) if ( numdeps .gt. kdim ) numdeps=kdim do 35 n=1,numdeps glub(n)=hall(n) cam(n)='XXXX' if ( glub(n) .lt. 10. ) then cam(n)='X ' elseif ( glub(n) .lt. 100. ) then cam(n)='XX ' elseif ( glub(n) .lt. 1000. ) then cam(n)='XXX ' endif inum=glub(n) call insertnum('X',inum,cam(n),4) 35 continue umask(1)='OCL_basin_AO.dat' umask(2)='OCL_basin_IO.dat' umask(3)='OCL_basin_PO.dat' umask(4)='OCL_basin_WO.dat' umask(5)='OCL_basin.dat' open(10,file='./basin_layer10ann_climatology.inf', * form='formatted', * status='old') read(10,'(a)')junk read(10,*) torm read(10,'(a)')junk read(10,*) woruw read(10,'(a)')junk read(10,*) aorv read(10,'(a)')junk read(10,'(a)')junk read(10,*) isetype read(10,'(a)')junk read(10,*) itime1 itime=itime1-1 read(10,'(a)')junk read(10,*) iwb1, iwb2 read(10,'(a)')junk read(10,'(a)') umask(5) read(10,'(a)')junk read(10,'(a)') bOCL read(10,'(a)')junk read(10,*) rsouth, rnorth read(10,'(a)')junk read(10,*) rw, re read(10,'(a)')junk read(10,*) (iup(i),ideep(i), i=1,10) read(10,*) (iup(i),ideep(i), i=11,20) read(10,*) (iup(i),ideep(i), i=21,30) read(10,'(a)')junk c read(10,*) ifst,ilst call clearstring(ctext,200) read(10,'(a)') ctext read(10,'(a)')junk read(10,*) pot read(10,'(a)')junk read(10,'(a)') infl read(10,'(a)')junk read(10,'(a6)') dotend read(10,'(a)')junk read(10,*) out_units read(10,'(a)')junk read(10,'(a)') out read(10,'(a)')junk read(10,'(a)') inse read(10,'(a)')junk read(10,'(a6)') doten1 read(10,'(a)')junk read(10,*) ifgrand read(10,'(a)')junk read(10,'(a)') ingrand read(10,'(a)')junk read(10,*) no_v read(10,'(a)') involkl read(10,'(a)')junk read(10,*) ivolobs read(10,'(a)') involyr close(10) c get first/last year and season/month call findlastchar(nsp,ctext,200) call spacefindl(nsp,icmax,ncomma,icomma,ctext) read(ctext(icomma(1,1):icomma(1,1)+icomma(1,2)-1),*) * ifst read(ctext(icomma(2,1):icomma(2,1)+icomma(2,2)-1),*) * ilst islst=itime1 if ( ncomma .ge. 3 ) then read(ctext(icomma(3,1):icomma(3,1)+icomma(3,2)-1),*) * islst write(6,*) 'islst',islst,ctext endif write(6,* )'ZZZ',ncomma,islst C How many layers were requested do i=1,30 if(iup(i).eq.0)then num_layers=i-1 goto 5555 endif enddo 5555 continue if(itime1.eq.6)then C monthly climatology baseyear=0. timestep=0. itime1=1 itime=0 else baseyear=1900.0 if(itime1.eq.5)then C mid of pentad timestep=2.5 islst=1 else C yearly, 6mons, 3mons, 1mons timestep=0.5/real(itime1) endif endif C Convert geographical coordinates into grid ones call geo2grid(rsouth,rdum,jsouth,idum,csouth,cdum) call geo2grid(rnorth,rdum,jnorth,idum,cnorth,cdum) C LAYER loop do 4000 iterla=1,num_layers nup=iup(iterla) ndeep=ideep(iterla) C middle depth of layer dila=0.5*(glub(nup)+glub(ndeep)) C open file with ODSQ GRAND MEAN anomalies of the layer if(ifgrand.eq.1)then cstr(1)=cam(nup) cstr(2)='-' cstr(3)=cam(ndeep) cstr(4)='m' call janameit(4,ingrand,'$',mess,cstr) call addCend(mess,200,0) call fileopen(mess,'rb+'//CHAR(0),10) call readlevel(1,1,1,1,1,idim,jdim,idim,jdim,idim,jdim, * kid,kid,10,f2) else call rinit2D(f2,idim,jdim,0.) endif C open file with ODSQ volumes of the layer if(no_v.ne.0)then cstr(1)=cam(nup) cstr(2)='-' cstr(3)=cam(ndeep) cstr(4)='m' call janameit(4,involkl,'$',mess,cstr) call addCend(mess,200,0) call fileopen(mess,'rb+'//CHAR(0),20) call readlevel(1,1,1,1,1,idim,jdim,idim,jdim,idim,jdim, * kid,kid,20,volk) else call rinit2D(volk,idim,jdim,0.) endif c BASIN loop do 3000 iocn=iwb1,iwb2 C Mask OCL basins to be excluded mess=umask(iocn) call basinOCLsubALL(fmask,idim,jdim,kdim,glub, * mess,bOCL,numdeps) c open output file for ITERLA-th layer cstr(1)=ox(iocn) cstr(2)=cam(nup) cstr(3)='-' cstr(4)=cam(ndeep) cstr(5)='m' cstr(6)=csouth cstr(7)=cnorth cstr(8)='.dat' call janameit(8,out,'$',mess,cstr) call addCend(mess,200,0) irod=18 open(irod,file=mess,status='unknown',form='formatted') c Initialize TIME loop if(itime.eq.4) then c pentads ilast=ilst-4 naTE=4 itime1=1 i6mod=0 else c yearly, 6mons, 3mons, monthly climatology ilast=ilst naTE=0 i6mod=itime endif c TIME loop do 2000 nfl=ifst, ilast call nameTIME(naTE,nfl,cyr,cyr5) c "subTIME" loop itimeloop=itime1 if ( nfl .eq. ilast ) itimeloop=islst write(6,*) 'timeloop',nfl,ilast,islst,itime1,itimeloop, * timestep do 2000 i6=1,itimeloop nz=0 volbzn=0. volobs=0. volprc=0. areabzn=0. addyear=(2*i6-1)*timestep outyear=baseyear+nfl+addyear revyear=int(1000.*addyear)+real(0.01*nfl) cstr(1)=cam(nup) cstr(2)='-' cstr(3)=cam(ndeep) cstr(4)='-' cstr(5)=cyr cstr(6)=cyr5 cstr(7)=c6mon(i6+i6mod) cstr(8)=dotend call janameit(8,infl,'$',mess,cstr) call addCend(mess,200,0) call fileopen(mess,'rb+'//CHAR(0),12) write(6,*) '6mon',inum,i6,i6mod,c6mon(i6+i6mod) if(isetype.gt.0)then cstr(8)=doten1 call janameit(8,inse,'$',mess,cstr) call addCend(mess,200,0) call fileopen(mess,'rb+'//CHAR(0),14) endif if(ivolobs.eq.1)then call janameit(7,involyr,'$',mess,cstr) call addCend(mess,200,0) call fileopen(mess,'rb+'//CHAR(0),21) endif call readlevel(1,1,1,1,1,idim,jdim,idim,jdim,idim,jdim, * kid,kid,12,f) if(isetype.gt.0)then call readlevel(1,1,1,1,1,idim,jdim,idim,jdim,idim,jdim, * kid,kid,14,f1) else call rinit2D(f1,idim,jdim,0.) endif if(ivolobs.eq.1)then call readlevel(1,1,1,1,1,idim,jdim,idim,jdim,idim,jdim, * kid,kid,21,fvol) else call rinit2D(fvol,idim,jdim,0.) endif c LATITUDE loop do 1900 j=jsouth,jnorth xlat=j-90.5 glat=xlat*pirat a1x1sq= a0 * odsqa(glat) if(rw.gt.re) re=hmsph+re irw=rw ire=re c do ro=rw,re do iro=irw,ire ro=iro ro=iro+0.5 ro1=trik(ro) ii=int(ro1+0.5) if(ii.gt.360) ii=ii-360 fval=f(ii,j) if(fval.gt.bmiss .and. fmask(ii,j,nup).gt.0.)then nz=nz+1 trim(nz)= (fval-f2(ii,j)) if(aorv.eq.1)then wtrim(nz)=a1x1sq elseif(aorv.eq.2)then wtrim(nz)=volk(ii,j) else wtrim(nz)=1. endif C 01/09/06 if there is no SE assign it to 0 if(f1(ii,j).gt.bmiss)then ferr(nz)=f1(ii,j) else ferr(nz)=0. endif volbzn=volbzn+volk(ii,j) areabzn=areabzn+a1x1sq volx=fvol(ii,j) if(volx.gt.bmiss) volobs=volobs+volx endif enddo 1900 continue !LATITUDE loop c compute and print c BASIN Total or Average for current basin and layer if (nz .gt. 0) then do i=1,nz y(i)=trim(i) yw(i)=wtrim(i) enddo C 02/02/06 this verson doesn't print VOLPRC, so I commented next line C if(volbzn.gt.0.) volprc=(volobs/volbzn)*100. if(torm.eq.1)then c basin total sum=0. do i=1,nz sum=sum+y(i) enddo sum = out_units * sum if(isetype.gt.0)then call basin_se(pot,nz,ferr(1),yw(1),1,isetype,x1,0.) x1 = out_units * x1 else x1=0. endif write(irod,'(f9.3,f8.1,2f15.6,i8,2f14.6,f8.2)') * outyear, * -1.*dila,sum,x1,nz,areabzn,volbzn, * revyear else c basin mean if(woruw.eq.1)then c weighted call trimmer(pot,nz,y(1),yw(1),3,x1,bmiss) x1 = out_units * x1 if(isetype.gt.0)then call basin_se(pot,nz,ferr(1),yw(1),3,isetype,temps,0.) temps = out_units * temps else temps = 0. endif write(irod,'(f9.3,f8.1,2f15.6,i8,2f14.6,f8.2)') * outyear, * -1.*dila,x1,temps,nz,areabzn,volbzn, * revyear else call trimmer(pot,nz,y(1),yw(1),2,x1,bmiss) x1 = out_units * x1 if(isetype.gt.0)then call basin_se(pot,nz,ferr(1),yw(1),2,isetype,temps,0.) temps= out_units * temps else temps = 0. endif write(irod,'(f9.3,f8.1,2f15.6,i8,2f14.6,f8.2)') * outyear, * -1.*dila,x1,temps,nz,areabzn,volbzn, * revyear endif endif endif nz=0 volobs=0. volbzn=0. volprc=0. areabzn=0. call rinit1D(trim,nrec,bmiss) call rinit1D(wtrim,nrec,bmiss) call rinit1D(ferr,nrec,0.) call onefileclose(12) if(isetype.gt.0) call onefileclose(14) if(ivolobs.eq.1) call onefileclose(21) 2000 continue !TIME loop close(irod) 3000 continue !BASIN loop if(ifgrand.eq.1) call onefileclose(10) if(no_v.ne.0) call onefileclose(20) 4000 continue !LAYER loop stop end c include 'antonov/sub_basinOCLrecent_tmod.f' c include 'antonov/sub_nameTIME.f' c include 'antonov/sub_fullname.f' c include 'antonov/sub_geo2OCLgrid.f' c include 'antonov/sub_trimmedmean.f' c include 'antonov/sub_init_array.f' c include 'antonov/sub_basin_se.f'