! implement an idea of Stuart's: compare (ux+vy)@iin with w<br>
! use doubled, subsampled grid to get divergence on w grid<br>
! properly accounts for partial bottom cells but not BBL (separate script)<br>
<br>
use subdomain_notides-mean.50s-20n.cdf<br>
use subdomain_notides-mean-w.50s-20n.cdf<br>
use acom3-cell-thickness.50s-20n.cdf		! from topog_acom3.1.dta.nc<br>
<br>
! make layer-integral transports:<br>
! first step just defines the first 36 values of e435 (u layer thickness in meters)<br>
let dhu = e435[d=acom3-cell-thickness.50s-20n, k=1:36]<br>
<br>
! second step makes transport in each layer (factor of 100 converts u from cm/s)<br>
! using these thicknesses takes the partial bottom cell into account<br>
let vdhu = v_30d[d=subdomain_notides-mean.50s-20n]*dhu/100<br>
let udhu = u_30d[d=subdomain_notides-mean.50s-20n]*dhu/100<br>
<br>
! double the grid for (u,v)<br>
define axis/x=90e:180:.25/unit=lon xhf<br>
define axis/y=50s:20n:`1/6`/unit=lat yhf<br>
define grid/like=v_30d[d=subdomain_notides-mean.50s-20n]/x=xhf/y=yhf ghf<br>
<br>
! transport components on the doubled grid<br>
let udhuhf = udhu[g=ghf]<br>
let vdhuhf = vdhu[g=ghf]<br>
<br>
! divergence doubled, then subsampled<br>
let divhf = udhuhf[x=@ddc] + vdhuhf[y=@ddc]<br>
let divhfrg = divhf[g=eta_30d[d=subdomain_notides-mean.50s-20n]]<br>
<br>
! integrate divergence down<br>
let divhfrgint = divhfrg[k=@rsum]<br>
<br>
! implicitly regrid div to w z-grid (thus lowering the value half a gridbox,<br>
! which is where the sum should accumulate to w).<br>
let diff = w_30d[d=subdomain_notides-mean-w.50s-20n]/100 - divhfrgint<br>