Changeset 4138764 in flexpart.git for src/netcdf_output_mod.f90

Timestamp:

Apr 7, 2021, 8:45:57 AM (3 years ago)

Author:

Sabine <sabine.eckhardt@…>

Branches:

dev

Children:

1228ef7

Parents:

03adec6 (diff), 759df5f (diff)
Note: this is a merge changeset, the changes displayed below correspond to the merge itself.
Use the (diff) links above to see all the changes relative to each parent.

Message:

Merge remote-tracking branch 'refs/remotes/origin/dev' into dev

File:

• src/netcdf_output_mod.f90 (modified) (8 diffs)

src/netcdf_output_mod.f90

@@ +1 @@
+! SPDX-FileCopyrightText: FLEXPART 1998-2019, see flexpart_license.txt
+! SPDX-License-Identifier: GPL-3.0-or-later
+
   !*****************************************************************************
   !                                                                            *
@@ -33 +36 @@
   use outg_mod,  only: outheight,oroout,densityoutgrid,factor3d,volume,&
                        wetgrid,wetgridsigma,drygrid,drygridsigma,grid,gridsigma,&
-                       area,arean,volumen, orooutn, areaeast, areanorth
+                       area,arean,volumen, orooutn, areaeast, areanorth, p0out, t0out
   use par_mod,   only: dep_prec, sp, dp, maxspec, maxreceptor, nclassunc,&
-                       unitoutrecept,unitoutreceptppt, nxmax,unittmp
+                       unitoutrecept,unitoutreceptppt, nxmax,unittmp, &
+                       write_p0t0
   use com_mod,   only: path,length,ldirect,bdate,ibdate,ibtime,iedate,ietime, &
                        loutstep,loutaver,loutsample,outlon0,outlat0,&
@@ -54 +58 @@
                        ioutputforeachrelease, iflux, mdomainfill, mquasilag, & 
                        nested_output, ipout, surf_only, linit_cond, &
-                       flexversion,mpi_mode,DRYBKDEP,WETBKDEP
+                       flexversion,mpi_mode,DRYBKDEP,WETBKDEP, ps, tt2
 
   use mean_mod
@@ -78 +82 @@
   integer, dimension(maxspec) :: specID,specIDppt, wdspecID,ddspecID
   integer, dimension(maxspec) :: specIDn,specIDnppt, wdspecIDn,ddspecIDn
+  integer :: psID, tt2ID
   integer                     :: timeID, timeIDn
   integer, dimension(6)       :: dimids, dimidsn
@@ -177 +182 @@
      call nf90_err(nf90_put_att(ncid, nf90_global, 'dyout', dyout))
   endif
-!       vertical levels stored in grid structure
+!   vertical levels stored in grid structure
 
   ! COMMAND file settings
@@ -385 +390 @@
   if (write_area) call nf90_err(nf90_def_var(ncid, 'area', nf90_float, &
        &(/ lonDimID, latDimID /), areaID))
+
+    ! surfarce pressure / temperature
+    if (write_p0t0) then
+       call nf90_err(nf90_def_var(ncid, 'pressure', nf90_float, &
+            &(/ lonDimID, latDimID, timeDimID /), psID))
+       call nf90_err(nf90_def_var(ncid, 'temperature', nf90_float, &
+            &(/ lonDimID, latDimID, timeDimID /), tt2ID))
+    end if
 
 
@@ -758 +771 @@
   ! real(sp)            :: wetgridtotal,wetgridsigmatotal
   ! real(sp)            :: drygridtotal,drygridsigmatotal
+    real :: ddx,ddy,p0h,t0h,p1,p2,p3,p4,rddx,rddy,xlon,xlat,ylat,xtn,ytn
+    integer :: i1,ixp,j1,jyp,j
 
   real, parameter     :: weightair=28.97
@@ -1038 +1053 @@
   end do
 
+
+    if (write_p0t0) then
+       ! Loop over all output grid cells
+       !********************************
+
+       do jjy=0,numygrid-1
+          do iix=0,numxgrid-1
+             p0h=0.
+             t0h=0.
+
+             ! Take 100 samples of the topography in every grid cell
+             !******************************************************
+
+             do j1=1,10
+                ylat=outlat0+(real(jjy)+real(j1)/10.-0.05)*dyout
+                yl=(ylat-ylat0)/dy
+                do i1=1,10
+                   xlon=outlon0+(real(iix)+real(i1)/10.-0.05)*dxout
+                   xl=(xlon-xlon0)/dx
+
+                   ! Determine the nest we are in
+                   !*****************************
+
+                   ngrid=0
+                   do j=numbnests,1,-1
+                      if ((xl.gt.xln(j)+eps).and.(xl.lt.xrn(j)-eps).and. &
+                           (yl.gt.yln(j)+eps).and.(yl.lt.yrn(j)-eps)) then
+                         ngrid=j
+                         goto 43
+                      endif
+                   end do
+43                 continue
+
+                   ! Determine (nested) grid coordinates and auxiliary parameters used for interpolation
+                   !*****************************************************************************
+
+                   if (ngrid.gt.0) then
+                      xtn=(xl-xln(ngrid))*xresoln(ngrid)
+                      ytn=(yl-yln(ngrid))*yresoln(ngrid)
+                      ix=int(xtn)
+                      jy=int(ytn)
+                      ddy=ytn-real(jy)
+                      ddx=xtn-real(ix)
+                   else
+                      ix=int(xl)
+                      jy=int(yl)
+                      ddy=yl-real(jy)
+                      ddx=xl-real(ix)
+                   endif
+                   ixp=ix+1
+                   jyp=jy+1
+                   rddx=1.-ddx
+                   rddy=1.-ddy
+                   p1=rddx*rddy
+                   p2=ddx*rddy
+                   p3=rddx*ddy
+                   p4=ddx*ddy
+
+                   p0h=p0h+p1*ps(ix ,jy,1,memind(1)) &
+                        + p2*ps(ixp,jy,1,memind(1)) &
+                        + p3*ps(ix ,jyp,1,memind(1)) &
+                        + p4*ps(ixp,jyp,1,memind(1))
+                   t0h=t0h+p1*tt2(ix ,jy,1,memind(1)) &
+                        + p2*tt2(ixp,jy,1,memind(1)) &
+                        + p3*tt2(ix ,jyp,1,memind(1)) &
+                        + p4*tt2(ixp,jyp,1,memind(1))
+                end do
+             end do
+
+             ! Divide by the number of samples taken
+             !**************************************
+             p0out(iix,jjy)=p0h/100.
+             t0out(iix,jjy)=t0h/100.
+          end do
+       end do
+
+       call nf90_err(nf90_put_var(ncid, psID, p0out, (/ 1,1,tpointer /), (/ numxgrid,numygrid,1 /)))
+       call nf90_err(nf90_put_var(ncid, tt2ID, t0out,(/ 1,1,tpointer /), (/ numxgrid,numygrid,1 /)))
+    end if
+
+
   ! Close netCDF file
   !**************************