[cdftools] 50/228: RD : add cdfpsi_level (streamfunction for each level) and cdfhdy (dynamic height computation)
Alastair McKinstry
mckinstry at moszumanska.debian.org
Fri Jun 12 08:21:28 UTC 2015
This is an automated email from the git hooks/post-receive script.
mckinstry pushed a commit to branch master
in repository cdftools.
commit 37018a0569cd721e649485edc3529a632538d0cd
Author: dussin <dussin at 1055176f-818a-41d9-83e1-73fbe5b947c5>
Date: Tue Jun 1 15:33:53 2010 +0000
RD : add cdfpsi_level (streamfunction for each level) and cdfhdy (dynamic height computation)
git-svn-id: http://servforge.legi.grenoble-inp.fr/svn/CDFTOOLS/trunk@326 1055176f-818a-41d9-83e1-73fbe5b947c5
---
Makefile | 9 +-
cdfhdy.f90 | 246 +++++++++++++++++++++++++++++++++++++++++++++++++++++++
cdfpsi_level.f90 | 166 +++++++++++++++++++++++++++++++++++++
3 files changed, 420 insertions(+), 1 deletion(-)
diff --git a/Makefile b/Makefile
index a85b2f6..b84b30d 100644
--- a/Makefile
+++ b/Makefile
@@ -31,7 +31,8 @@ EXEC = cdfmoy cdfmoyt cdfmoy_sp cdfstd cdfmoy_sal2_temp2 cdfmoy_annual cdfmoy_
cdfprofile cdfwhereij cdffindij cdfweight cdfmaxmoc cdfcensus cdfzoom cdfmax cdfmax_sp cdfprobe \
bimgmoy4 bimgcaltrans cdf16bit cdfvita cdfconvert cdfflxconv cdfclip cdfsstconv cdfstrconv cdfbathy cdfvar cdfmkmask-zone\
cdfcsp cdfcoloc cdfmltmask cdfstatcoord cdfpolymask cdfsmooth cdfmkmask cdfdifmask\
- cdfkempemekeepe cdfbci cdfbti cdfnrjcomp cdfcofdis cdfsections cdfnorth_unfold cdfovide cdfmppini
+ cdfkempemekeepe cdfbci cdfbti cdfnrjcomp cdfcofdis cdfsections cdfnorth_unfold cdfovide cdfmppini\
+ cdfpsi_level cdfhdy
all: $(EXEC)
@@ -202,6 +203,9 @@ cdfbti: cdfio.o cdfbti.f90
cdfnrjcomp: cdfio.o cdfnrjcomp.f90
$(F90) cdfnrjcomp.f90 -o cdfnrjcomp cdfio.o $(FFLAGS)
+cdfhdy: cdfio.o eos.o cdfhdy.f90
+ $(F90) cdfhdy.f90 -o cdfhdy cdfio.o eos.o $(FFLAGS)
+
## Transport programs
cdfmhst: cdfio.o cdfmhst.f90
$(F90) cdfmhst.f90 -o cdfmhst cdfio.o $(FFLAGS)
@@ -239,6 +243,9 @@ cdfpsi-austral: cdfio.o cdfpsi-austral.f90
cdfpsi-austral-ssh: cdfio.o cdfpsi-austral-ssh.f90
$(F90) cdfpsi-austral-ssh.f90 -o cdfpsi-austral-ssh cdfio.o $(FFLAGS)
+cdfpsi_level: cdfio.o cdfpsi_level.f90
+ $(F90) cdfpsi_level.f90 -o cdfpsi_level cdfio.o $(FFLAGS)
+
cdftransportiz: cdfio.o cdftransportiz.f90
$(F90) cdftransportiz.f90 -o cdftransportiz cdfio.o $(FFLAGS)
diff --git a/cdfhdy.f90 b/cdfhdy.f90
new file mode 100644
index 0000000..3e37c71
--- /dev/null
+++ b/cdfhdy.f90
@@ -0,0 +1,246 @@
+PROGRAM cdfhdy
+ !!-------------------------------------------------------------------
+ !! *** PROGRAM cdfhdy ***
+ !!
+ !! ** Purpose: Compute dynamical height anomaly field from gridT file
+ !! Store the results on a 2D cdf file.
+ !!
+ !!
+ !! history:
+ !! Original : J.M. Molines (Nov 2004 ) for ORCA025
+ !! J.M. Molines Apr 2005 : use modules
+ !!-------------------------------------------------------------------
+ !! $Rev: 256 $
+ !! $Date: 2009-07-21 17:49:27 +0200 (mar. 21 juil. 2009) $
+ !! $Id: cdfsig0.f90 256 2009-07-21 15:49:27Z molines $
+ !!--------------------------------------------------------------
+ !! * Modules used
+ USE cdfio
+ USE eos
+
+ !! * Local variables
+ IMPLICIT NONE
+ INTEGER :: jk, jt !: dummy loop index
+ INTEGER :: ierr !: working integer
+ INTEGER :: narg, iargc !:
+ INTEGER :: npiglo,npjglo, npk, npt !: size of the domain
+ INTEGER, DIMENSION(1) :: ipk, & !: outptut variables : number of levels,
+ & id_varout !: ncdf varid's
+ real(KIND=4) , DIMENSION (:,:), ALLOCATABLE :: ztemp, zsal ,& !: Array to read a layer of data
+ & ztemp0, zsal0 ,& !: reference density
+ & zsig0 , & !: potential density (sig-0)
+ & zsig , & !: potential density (sig-0)
+ & zmask , & !: 2D mask at current level
+ & zhdy, zterm, zdep, zdepth, zssh
+ REAL(KIND=4),DIMENSION(:),ALLOCATABLE :: tim
+
+ CHARACTER(LEN=256) :: cfilet , cdum, cfileout='cdfhdy.nc', cmask='mask.nc' !:
+ CHARACTER(LEN=256) :: coordzgr='mesh_zgr.nc'
+
+ TYPE(variable) , DIMENSION(1) :: typvar !: structure for attributes
+
+ INTEGER :: ncout
+ INTEGER :: istatus
+ INTEGER :: zlev1, zlev2
+ INTEGER, DIMENSION (2) :: ismin, ismax
+ REAL(KIND=4) :: sigmin, sigmax, rau0=1000.
+
+
+ !! Read command line
+ narg= iargc()
+ IF ( narg .LT. 3 ) THEN
+ PRINT *,' Usage : cdfhdy gridT level1 level2 '
+ PRINT *,' integrates from level1 (usually surface) to level2, level2 greater than level1 '
+ PRINT *,' reference is the sea surface, mask.nc and mesh_zgr.nc must be in your directory'
+ PRINT *,' Output on cdfhdy.nc, variable sohdy'
+ STOP
+ ENDIF
+
+ CALL getarg (1, cfilet)
+ CALL getarg (2, cdum) ; READ(cdum,*) zlev1
+ CALL getarg (3, cdum) ; READ(cdum,*) zlev2
+ npiglo= getdim (cfilet,'x')
+ npjglo= getdim (cfilet,'y')
+ npk = getdim (cfilet,'depth')
+ npt = getdim (cfilet,'time')
+
+ ipk(:)= 1
+ typvar(1)%name= 'sohdy'
+ typvar(1)%units='m'
+ typvar(1)%missing_value=0.
+ typvar(1)%valid_min= -100.
+ typvar(1)%valid_max= 100.
+ typvar(1)%long_name='Dynamical height anomaly'
+ typvar(1)%short_name='sohdy'
+ typvar(1)%online_operation='N/A'
+ typvar(1)%axis='TZYX'
+
+
+ PRINT *, 'npiglo=', npiglo
+ PRINT *, 'npjglo=', npjglo
+ PRINT *, 'npk =', npk
+ PRINT *, 'npt =', npt
+
+ ALLOCATE (ztemp0(npiglo,npjglo), zsal0(npiglo,npjglo), zsig0(npiglo,npjglo) ,zmask(npiglo,npjglo))
+ ALLOCATE (ztemp(npiglo,npjglo), zsal(npiglo,npjglo), zsig(npiglo,npjglo) , zhdy(npiglo,npjglo), zterm(npiglo,npjglo))
+ ALLOCATE (zdep(npiglo,npjglo), zdepth(npiglo,npjglo), zssh(npiglo,npjglo))
+ ALLOCATE (tim(npt))
+
+ ! create output fileset
+
+ ncout =create(cfileout, cfilet, npiglo,npjglo,npk)
+
+ ierr= createvar(ncout ,typvar,1, ipk,id_varout )
+ ierr= putheadervar(ncout, cfilet,npiglo, npjglo,npk)
+ tim=getvar1d(cfilet,'time_counter',npt)
+ ierr=putvar1d(ncout,tim,npt,'T')
+
+ ! Temperature and salinity for reference profile
+ ztemp0(:,:)=0.
+ zsal0(:,:)=35.
+
+ zmask(:,:) = getvar(cmask, 'tmask', zlev2, npiglo, npjglo)
+
+ DO jt=1,npt
+ PRINT *,' TIME = ', jt, tim(jt)/86400.,' days'
+
+ zhdy(:,:) = 0.
+ zdepth(:,:) = 0.
+
+ DO jk = zlev1, zlev2
+
+ zdep(:,:) = getvar(coordzgr, 'e3t_ps', jk,npiglo,npjglo,ldiom=.true.)
+
+ ! total depth at current level (used for computation of rho in situ)
+ zdepth(:,:) = zdepth(:,:) + zdep(:,:)
+
+ ztemp(:,:)= getvar(cfilet, 'votemper', jk ,npiglo, npjglo,ktime=jt)
+ zsal(:,:) = getvar(cfilet, 'vosaline', jk ,npiglo, npjglo,ktime=jt)
+
+ CALL eos_insitu( ztemp0, zsal0, zdepth, npiglo, npjglo, zsig0 )
+ CALL eos_insitu( ztemp, zsal, zdepth, npiglo, npjglo, zsig )
+
+ PRINT *, 'max of ref profile for level', jk ,'is ', MAXVAL(zsig0)
+
+ ! we compute the term of the integral : (1/g) *10e4 * sum [ delta * dz ]
+ ! with delta = (1/rho - 1/rho0)
+ ! 10e4 factor is conversion decibar/pascal
+ !
+ zterm = ( ( 1. / ( rau0 + zsig(:,:) ) ) - ( 1. / ( rau0 + zsig0(:,:) ) ) ) * 10000. * zdep / 9.81
+ ! in land, it seems appropriate to stop the computation
+ WHERE(zsal == 0 ) zterm = 0
+
+ zhdy(:,:) = zhdy(:,:) + zterm(:,:)
+
+ END DO ! loop to next level
+
+ ! we mask with the last level of the integral
+ zhdy(:,:) = zhdy(:,:) * zmask(:,:)
+
+ ierr = putvar(ncout, id_varout(1) ,zhdy, 1,npiglo, npjglo,ktime=jt)
+
+ END DO ! next time frame
+
+ istatus = closeout(ncout)
+
+CONTAINS
+
+SUBROUTINE eos_insitu( ptem, psal, pdepth, jpiglo, jpjglo, prd )
+ !!----------------------------------------------------------------------
+ !! *** ROUTINE eos_insitu ***
+ !!
+ !! ** Purpose : Compute the in situ density (ratio rho/rau0) from
+ !! potential temperature and salinity using an equation of state
+ !! defined through the namelist parameter nn_eos.
+ !!
+ !! ** Method :
+ !! nn_eos = 0 : Jackett and McDougall (1994) equation of state.
+ !! the in situ density is computed directly as a function of
+ !! potential temperature relative to the surface (the opa t
+ !! variable), salt and pressure (assuming no pressure variation
+ !! along geopotential surfaces, i.e. the pressure p in decibars
+ !! is approximated by the depth in meters.
+ !! prd(t,s,p) = ( rho(t,s,p) - rau0 ) / rau0
+ !! with pressure p decibars
+ !! potential temperature t deg celsius
+ !! salinity s psu
+ !! reference volumic mass rau0 kg/m**3
+ !! in situ volumic mass rho kg/m**3
+ !! in situ density anomalie prd no units
+ !! Check value: rho = 1060.93298 kg/m**3 for p=10000 dbar,
+ !! t = 40 deg celcius, s=40 psu
+ !! prd(t,s) = rn_beta * s - rn_alpha * tn - 1.
+ !! Note that no boundary condition problem occurs in this routine
+ !! as (ptem,psal) are defined over the whole domain.
+ !!
+ !! ** Action : compute prd , the in situ density (no units)
+ !!
+ !! References : Jackett and McDougall, J. Atmos. Ocean. Tech., 1994
+ !!----------------------------------------------------------------------
+ INTEGER, INTENT(in ) :: jpiglo, jpjglo
+ REAL(4), DIMENSION(jpiglo,jpjglo), INTENT(in ) :: ptem ! potential temperature [Celcius]
+ REAL(4), DIMENSION(jpiglo,jpjglo), INTENT(in ) :: psal ! salinity [psu]
+ REAL(4), DIMENSION(jpiglo,jpjglo), INTENT(in ) :: pdepth ! depth [m]
+ REAL(4), DIMENSION(jpiglo,jpjglo), INTENT( out) :: prd ! in situ density
+ !!
+ INTEGER :: ji, jj, jk ! dummy loop indices
+ INTEGER :: jpkm1
+ REAL(4) :: zt , zs , zh , zsr ! temporary scalars
+ REAL(4) :: zr1, zr2, zr3, zr4 ! - -
+ REAL(4) :: zrhop, ze, zbw, zb ! - -
+ REAL(4) :: zd , zc , zaw, za ! - -
+ REAL(4) :: zb1, za1, zkw, zk0 ! - -
+ REAL(4) :: zrau0r ! - -
+ REAL(4), DIMENSION(jpiglo,jpjglo) :: zws ! temporary workspace
+ INTEGER :: nn_eos = 0 !: = 0/1/2 type of eq. of state and Brunt-Vaisala frequ.
+ REAL(4) :: rn_alpha = 2.0e-4 !: thermal expension coeff. (linear equation of state)
+ REAL(4) :: rn_beta = 7.7e-4 !: saline expension coeff. (linear equation of state)
+
+ REAL(4) :: ralpbet !: alpha / beta ratio
+ !!----------------------------------------------------------------------
+
+ zrau0r = 1.e0 / rau0
+ zws(:,:) = SQRT( ABS( psal(:,:) ) )
+ !
+ DO jj = 1, jpjglo
+ DO ji = 1, jpiglo
+ zt = ptem (ji,jj)
+ zs = psal (ji,jj)
+ zh = pdepth(ji,jj) ! depth
+ zsr= zws (ji,jj) ! square root salinity
+ !
+ ! compute volumic mass pure water at atm pressure
+ zr1= ( ( ( ( 6.536332e-9*zt-1.120083e-6 )*zt+1.001685e-4)*zt &
+ & -9.095290e-3 )*zt+6.793952e-2 )*zt+999.842594
+ ! seawater volumic mass atm pressure
+ zr2= ( ( ( 5.3875e-9*zt-8.2467e-7 ) *zt+7.6438e-5 ) *zt &
+ & -4.0899e-3 ) *zt+0.824493
+ zr3= ( -1.6546e-6*zt+1.0227e-4 ) *zt-5.72466e-3
+ zr4= 4.8314e-4
+ !
+ ! potential volumic mass (reference to the surface)
+ zrhop= ( zr4*zs + zr3*zsr + zr2 ) *zs + zr1
+ !
+ ! add the compression terms
+ ze = ( -3.508914e-8*zt-1.248266e-8 ) *zt-2.595994e-6
+ zbw= ( 1.296821e-6*zt-5.782165e-9 ) *zt+1.045941e-4
+ zb = zbw + ze * zs
+ !
+ zd = -2.042967e-2
+ zc = (-7.267926e-5*zt+2.598241e-3 ) *zt+0.1571896
+ zaw= ( ( 5.939910e-6*zt+2.512549e-3 ) *zt-0.1028859 ) *zt - 4.721788
+ za = ( zd*zsr + zc ) *zs + zaw
+ !
+ zb1= (-0.1909078*zt+7.390729 ) *zt-55.87545
+ za1= ( ( 2.326469e-3*zt+1.553190)*zt-65.00517 ) *zt+1044.077
+ zkw= ( ( (-1.361629e-4*zt-1.852732e-2 ) *zt-30.41638 ) *zt + 2098.925 ) *zt+190925.6
+ zk0= ( zb1*zsr + za1 )*zs + zkw
+ !
+ ! masked in situ density anomaly
+ prd(ji,jj) = ( zrhop / ( 1.0 - zh / ( zk0 - zh * ( za - zh * zb ) ) ) &
+ & - rau0 ) ! * zrau0r ! * tmask(ji,jj)
+ END DO
+ END DO
+END SUBROUTINE eos_insitu
+
+END PROGRAM cdfhdy
diff --git a/cdfpsi_level.f90 b/cdfpsi_level.f90
new file mode 100644
index 0000000..10a882a
--- /dev/null
+++ b/cdfpsi_level.f90
@@ -0,0 +1,166 @@
+PROGRAM cdfpsi_level
+ !!-------------------------------------------------------------------
+ !! *** PROGRAM cdfpsi_level ***
+ !!
+ !! ** Purpose : Compute Stream Function for each level
+ !! PARTIAL STEPS
+ !!
+ !! ** Method : Compute the 3D fields ztrpu, ztrpv
+ !! as the integral on the vertical of u, v on their
+ !! respective points.
+ !! Then integrate from south to north : ==> psiu
+ !! Then integrate from West to East : ==> psiv
+ !! (should be almost the same (if no error ))
+ !! Default (appropriate for global model): output psiu;
+ !! normalizes the values setting psi (jpi,jpj) = 0
+ !! If option "V" is given as last argument, output psiv,
+ !! normalizes values setting psi(jpi,1) = 0.
+ !! This is appropriate for North Atlantic
+ !!
+ !! history ;
+ !! Original : J.M. Molines (May 2005 )
+ !!-------------------------------------------------------------------
+ !! $Rev: 256 $
+ !! $Date: 2009-07-21 17:49:27 +0200 (mar. 21 juil. 2009) $
+ !! $Id: cdfpsi.f90 256 2009-07-21 15:49:27Z molines $
+ !!--------------------------------------------------------------
+ !! * Modules used
+ USE cdfio
+
+ !! * Local variables
+ IMPLICIT NONE
+ INTEGER :: ji,jj,jk !: dummy loop index
+ INTEGER :: ierr !: working integer
+ INTEGER :: narg, iargc !: command line
+ INTEGER :: npiglo,npjglo, npk !: size of the domain
+ INTEGER :: ncout
+ INTEGER, DIMENSION(1) :: ipk, id_varout !
+
+ REAL(KIND=4), DIMENSION (:,:), ALLOCATABLE :: zmask, e1v, e3v , zv !: mask, metrics
+ REAL(KIND=4), DIMENSION (:,:), ALLOCATABLE :: e2u, e3u , zu !: mask, metrics
+ REAL(KIND=4), DIMENSION (:,:), ALLOCATABLE :: glamf, gphif
+ REAL(KIND=4) ,DIMENSION(1) :: tim
+
+ REAL(KIND=8), DIMENSION (:,:), ALLOCATABLE :: ztrpu, ztrpv, psiu, psiv
+
+ CHARACTER(LEN=256) :: cfileu ,cfilev, cfileoutnc='psi_level.nc'
+ CHARACTER(LEN=256) :: coordhgr='mesh_hgr.nc', coordzgr='mesh_zgr.nc', cmask='mask.nc'
+ CHARACTER(LEN=1) :: coption
+ CHARACTER(LEN=256) :: cdep
+
+ TYPE(variable), DIMENSION(1) :: typvar !: structure for attributes
+
+ INTEGER :: istatus
+
+ ! constants
+
+ !! Read command line and output usage message if not compliant.
+ narg= iargc()
+ IF ( narg == 0 ) THEN
+ PRINT *,' Usage : cdfpsi_level Ufile Vfile <V> (optional argument)'
+ PRINT *,' Computes the barotropic stream function as the integral of the transport'
+ PRINT *,' PARTIAL CELLS VERSION'
+ PRINT *,' Files mesh_hgr.nc, mesh_zgr.nc ,mask.nc must be in te current directory'
+ PRINT *,' Output on psi_level.nc, variables sobarstf on f-points'
+ PRINT *,' Default works well for a global ORCA grid. use V 3rdargument for North Atlantic'
+ STOP
+ ENDIF
+
+ CALL getarg (1, cfileu )
+ CALL getarg (2, cfilev )
+ CALL getarg (3, coption )
+
+ npiglo= getdim (cfileu,'x')
+ npjglo= getdim (cfileu,'y')
+ npk = getdim (cfileu,'depth')
+
+ ! define new variables for output ( must update att.txt)
+ typvar(1)%name= 'sobarstf'
+ typvar(1)%units='m3/s'
+ typvar(1)%missing_value=0.
+ typvar(1)%valid_min= -300.e6
+ typvar(1)%valid_max= 300.e6
+ typvar(1)%long_name='Barotropic_Stream_Function'
+ typvar(1)%short_name='sobarstf'
+ typvar(1)%online_operation='N/A'
+ typvar(1)%axis='TZYX'
+ ipk(1) = npk ! 3D ( X, Y , Z, T )
+
+ PRINT *, 'npiglo=', npiglo
+ PRINT *, 'npjglo=', npjglo
+ PRINT *, 'npk =', npk
+ IF ( coption == 'V') PRINT *, ' Use psiv (ex. North Atlantic case)'
+
+ ! Allocate arrays
+ ALLOCATE ( zmask(npiglo,npjglo) )
+ ALLOCATE ( e1v(npiglo,npjglo),e3v(npiglo,npjglo))
+ ALLOCATE ( e2u(npiglo,npjglo),e3u(npiglo,npjglo))
+ ALLOCATE ( zu(npiglo,npjglo),ztrpu(npiglo,npjglo), psiu(npiglo,npjglo) )
+ ALLOCATE ( zv(npiglo,npjglo),ztrpv(npiglo,npjglo), psiv(npiglo,npjglo))
+ ALLOCATE ( glamf(npiglo,npjglo), gphif(npiglo,npjglo))
+
+ glamf(:,:) = getvar(coordhgr, 'glamf',1,npiglo,npjglo)
+ gphif(:,:) = getvar(coordhgr, 'gphif',1,npiglo,npjglo)
+
+ ! create output fileset
+ ncout =create(cfileoutnc, cfileu, npiglo,npjglo,npk)
+ ierr= createvar(ncout ,typvar,1, ipk,id_varout )
+ ierr= putheadervar(ncout , cfileu, npiglo, npjglo, npk)
+! ierr= putheadervar(ncout , cfileu, npiglo, npjglo, npk,cdep=cdep)
+! ierr= putheadervar(ncout, cfileu,npiglo, npjglo,1,glamf, gphif)
+ tim=getvar1d(cfileu,'time_counter',1)
+ ierr=putvar1d(ncout,tim,1,'T')
+
+ e1v(:,:) = getvar(coordhgr, 'e1v', 1,npiglo,npjglo)
+ e2u(:,:) = getvar(coordhgr, 'e2u', 1,npiglo,npjglo)
+
+ ztrpu(:,:)= 0.d0
+ ztrpv(:,:)= 0.d0
+
+ DO jk = 1,npk
+
+ zmask(:,:) = getvar(cmask, 'fmask', jk,npiglo,npjglo)
+ ! get rid of the free-slip/no-slip condition
+ WHERE ( zmask >= 2 ) zmask = 1
+
+ PRINT *,'level ',jk
+ IF ( coption == 'V' ) THEN
+ zv(:,:)= getvar(cfilev, 'vomecrty', jk ,npiglo,npjglo)
+ e3v(:,:) = getvar(coordzgr, 'e3v_ps', jk,npiglo,npjglo, ldiom=.true.)
+ ztrpv(:,:) = zv(:,:)*e1v(:,:)*e3v(:,:) ! meridional transport of each grid cell
+ ELSE
+ ! Get zonal velocity at jk
+ zu(:,:)= getvar(cfileu, 'vozocrtx', jk ,npiglo,npjglo)
+ ! get e3v at level jk
+ e3u(:,:) = getvar(coordzgr, 'e3u_ps', jk,npiglo,npjglo, ldiom=.true.)
+ ! integrates vertically
+ ztrpu(:,:) = zu(:,:)*e2u(:,:)*e3u(:,:) ! zonal transport of each grid cell
+ ENDIF
+
+ IF (coption == 'V' ) THEN
+ ! integrate zonally from east to west
+ psiv(npiglo,:)= 0.0
+ DO ji=npiglo-1,1,-1
+ psiv(ji,:) = psiv(ji+1,:) - ztrpv(ji,:) ! psi at f point
+ END DO
+ psiv(:,:) = psiv(:,:) *zmask(:,:)
+ !ierr = putvar(ncout, id_varout(1) ,REAL(psiv), jk, npiglo, npjglo)
+ ierr = putvar(ncout, id_varout(1) ,REAL(ztrpv), jk, npiglo, npjglo)
+
+ ELSE
+ ! integrate from the south to the north with zonal transport
+ psiu(:,:) = 0.d0
+
+ DO jj = 2, npjglo
+ psiu(:,jj) = psiu(:,jj-1) - ztrpu(:,jj) ! psi at f point
+ END DO
+ psiu(:,:) = (psiu(:,:) -psiu(npiglo,npjglo) ) * zmask(:,:)
+ ierr = putvar(ncout, id_varout(1) ,REAL(psiu), jk, npiglo, npjglo)
+ !ierr = putvar(ncout, id_varout(1) ,REAL(ztrpu), jk, npiglo, npjglo)
+ ENDIF
+
+ END DO ! loop to next level
+
+ istatus = closeout (ncout)
+
+ END PROGRAM cdfpsi_level
--
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-science/packages/cdftools.git
More information about the debian-science-commits
mailing list