subroutine set_depths c c======================================================================= c === c This routine returns the depth at each grid point along a slab === c of constant J in the hybrid coordinate system. === c === c ------ === c Input: === c ------ === c === c Common Blocks: === c === c /HYBRID/ === c === c HAVE Total thickness of reference === C depths in lower levels. (real; cm) === c HC Total thickness of reference === c depths in upper levels. (real; cm) === c KC Coordinate interface level (integer) === c REFZ Reference depths (real array; cm) === #if defined dblsigma c ZC1 Minimum depth of the coordinate === c interface. (real; cm) === c ZC2 Maximum depth of the coordinate === c interface. (real; cm) === c ZREF Coordinate interface ref. depth. (real; cm) === c ZSLOPE Coordinate interface slope param. (real; cm) === #endif c === c /VERTSLAB/ === c === c HDV The bottom depth at velocity === c grid points. (real array; cm) === c === c ------- === c Output: === c ------- === c === c Common Blocks: === c === c /VERTSLAB/ === c === c HD Bottom depth at tracer grid points. (real array; cm) === c HVAV Reciprocal thickness, velocity points. (real array; cm) === c === c /VERTICAL/ === c === c HVZ Depths at uv 3D grid points. (real array; cm) === c HTZ Depths at tracer 3D grid points. (real array; cm) === c === c Calls: GET_THICK === c === c======================================================================= c c----------------------------------------------------------------------- c Define global data. c----------------------------------------------------------------------- c #include #include #include #include #include #include c c----------------------------------------------------------------------- c Define local data. c----------------------------------------------------------------------- c logical first integer i,j,k FLOAT * thickv(imt,km,0:1),dum(imt,kmp1) #ifdef dblsigma FLOAT * fac,favg,fdif,fv,hvmfv,rhave,rhc save first,rhave,rhc,fac,favg,fdif #else FLOAT * rhave,sigma save first,rhave #endif data first /.true./ c c======================================================================= c Begin executable code. c======================================================================= c #ifdef dblsigma c c----------------------------------------------------------------------- c Set often used interface parameters, on first call. c----------------------------------------------------------------------- c if(first) then favg=p5*(zc1+zc2) fdif=p5*(zc2-zc1) if(fdif.ne.c0) then fac=zslope/fdif else fac=c0 endif if (kc.gt.0) then rhc=c1/hc else rhc=c0 endif rhave=c1/have first=.false. endif c c----------------------------------------------------------------------- c Compute depths at each horizontal grid point. c----------------------------------------------------------------------- c do 30 j=1,jmt do 30 i=1,imt fv=favg+fdif*tanh(fac*(hdv(i,j)-zref)) hvmfv=hdv(i,j)-fv c c----------------------------------------------------------------------- c In upper sigma system. c----------------------------------------------------------------------- c do 10 k=1,kc hvz(i,k,j)=refz(k)*rhc*fv 10 continue c c----------------------------------------------------------------------- c In lower sigma system. c----------------------------------------------------------------------- c do 20 k=kc+1,km hvz(i,k,j)=fv+(refz(k)-hc)*rhave*hvmfv 20 continue 30 continue #else c c----------------------------------------------------------------------- c On first call, set reciprocal of have. c----------------------------------------------------------------------- c if(first) then rhave=c1/have first=.false. endif c c----------------------------------------------------------------------- c Determine depths in flat levels. c----------------------------------------------------------------------- c do 40 j=1,jmt do 20 k=1,kc do 10 i=1,imt hvz(i,k,j)=refz(k) 10 continue 20 continue c c----------------------------------------------------------------------- c Determine depths in "sigma" levels. c----------------------------------------------------------------------- c do 40 k=kc+1,km sigma=(refz(k)-hc)*rhave do 30 i=1,imt hvz(i,k,j)=hc+sigma*(hdv(i,j)-hc) 30 continue 40 continue #endif c c----------------------------------------------------------------------- c Determine depths at tracer points. c Depths at center of tracer boxes htz and bottom depth hd c----------------------------------------------------------------------- c do 90 j=2,jmt call get_thick(hvz(1,1,j-1),thickv(1,1,0),dum,kmp1) call get_thick(hvz(1,1, j),thickv(1,1,1),dum,kmp1) do 80 i=2,imt do 50 k=1,km htz(i,k,j)=p25*(thickv(i-1,k,0)+thickv( i,k,0)+ & thickv(i-1,k,1)+thickv( i,k,1)) 50 continue dum(i,1)=p5*htz(i,1,j) #ifndef barotropic do 60 k=2,km dum(i,k)=p5*(htz(i,k-1,j)+htz(i,k,j))+dum(i,k-1) 60 continue #endif hd(i,j)=dum(i,km)+p5*htz(i,km,j) do 70 k=1,km htz(i,k,j)=dum(i,k) 70 continue 80 continue 90 continue c c----------------------------------------------------------------------- c Determine depths at tracer edges. j=1,i=1 c Depths at center of tracer boxes htz and bottom depth hd c----------------------------------------------------------------------- c #ifdef cyclic do 100 i=2,imt hd(i,1)=hd(i,2) do 100 k=1,km htz(i,k,1)=htz(i,k,2) 100 continue do 110 j=1,jmt hd(1,j)=hd(imtm2,j) hd(imt,j)=hd(2,j) do 110 k=1,km htz(1,k,j)=htz(imtm2,k,j) htz(imt,k,j)=htz(2,k,j) 110 continue c #else do 100 j=2,jmt hd(1,j)=hd(2,j) do 100 k=1,km htz(1,k,j)=htz(2,k,j) 100 continue do 110 i=2,imt hd(i,1)=hd(i,2) do 110 k=1,km htz(i,k,1)=htz(i,k,2) 110 continue hd(1,1)=p5*(hd(1,2)+hd(2,1)) do 120 k=1,km htz(1,k,1)=p5*(htz(1,k,2)+htz(2,k,1)) 120 continue c #endif c----------------------------------------------------------------------- c Compute total thickness of the water column over velocity grid. c----------------------------------------------------------------------- c do 130 j = 1, jmt do 130 i = 1, imt hvav(i,j) = c1/hdv(i,j) 130 continue c return end