Changeset 37 for branches/petra/src/interpol_rain_nests.f90

Timestamp:

Apr 22, 2015, 3:42:35 PM (9 years ago)

Author:

pesei

Message:

Wet dep quick fix and other small changes. Wet depo quick fix not final yet.

File:

• branches/petra/src/interpol_rain_nests.f90 (modified) (7 diffs)

branches/petra/src/interpol_rain_nests.f90

@@ -1 +1 @@
 !**********************************************************************
-! Copyright 1998,1999,2000,2001,2002,2005,2007,2008,2009,2010         *
+! Copyright 1998-2015                                                 *
 ! Andreas Stohl, Petra Seibert, A. Frank, Gerhard Wotawa,             *
 ! Caroline Forster, Sabine Eckhardt, John Burkhart, Harald Sodemann   *
@@ -21 +21 @@
 
 subroutine interpol_rain_nests(yy1,yy2,yy3,nxmaxn,nymaxn,nzmax, &
-       maxnests,ngrid,nxn,nyn,memind,xt,yt,level,itime1,itime2,itime, &
-       yint1,yint2,yint3)
+       maxnests,ngrid,nxn,nyn, &
+       memind,xt,yt,level,itime1,itime2,itime,yint1,yint2,yint3, &
+       icbot,ictop)
   !                                i   i   i    i      i      i
   !   i       i    i   i    i    i  i    i     i      i      i
@@ -44 +45 @@
   !                                                                           *
   !     15 March 2000                                                         *
+  !                                                                           *
+  !****************************************************************************
+  ! Petra Seibert, 2011/2012-2015:
+  !  Add interpolation of cloud bottom and cloud thickness
+  !  for fix to SE's new wet scavenging scheme
   !                                                                           *
   !****************************************************************************
@@ -71 +77 @@
   !****************************************************************************
 
+  use com_mod, only: icloudbotn, icloudthckn
+  use par_mod, only: icmv
+
   implicit none
 
   integer :: maxnests,ngrid
-  integer :: nxn(maxnests),nyn(maxnests),nxmaxn,nymaxn,nzmax,memind(2)
-  integer :: m,ix,jy,ixp,jyp,itime,itime1,itime2,level,indexh
+  integer :: nxn(maxnests),nyn(maxnests),nxmaxn,nymaxn,nzmax,memind(2),m
+  integer :: ix,jy,ixp,jyp
+  integer :: itime,itime1,itime2,level,indexh,ip1,ip2,ip3,ip4
+  integer :: icbot,ictop,icthck,ipsum
   real :: yy1(0:nxmaxn-1,0:nymaxn-1,nzmax,2,maxnests)
   real :: yy2(0:nxmaxn-1,0:nymaxn-1,nzmax,2,maxnests)
   real :: yy3(0:nxmaxn-1,0:nymaxn-1,nzmax,2,maxnests)
-  real :: ddx,ddy,rddx,rddy,dt1,dt2,dt,y1(2),y2(2),y3(2)
+  real :: ddx,ddy,rddx,rddy,dt1,dt2,dt,y1(2),y2(2),y3(2),cbot(2),cthck(2)
   real :: xt,yt,yint1,yint2,yint3,p1,p2,p3,p4
-
-
 
   ! If point at border of grid -> small displacement into grid
   !***********************************************************
 
-  if (xt.ge.(real(nxn(ngrid)-1)-0.0001)) &
-       xt=real(nxn(ngrid)-1)-0.0001
-  if (yt.ge.(real(nyn(ngrid)-1)-0.0001)) &
-       yt=real(nyn(ngrid)-1)-0.0001
+  if (xt.ge.(real(nxn(ngrid)-1)-0.0001)) xt=real(nxn(ngrid)-1)-0.0001
+  if (yt.ge.(real(nyn(ngrid)-1)-0.0001)) yt=real(nyn(ngrid)-1)-0.0001
 
 
@@ -119 +126 @@
   !***********************
 
+  memind_loop: &
   do m=1,2
+
     indexh=memind(m)
 
@@ -134 +143 @@
          + p3*yy3(ix ,jyp,level,indexh,ngrid) &
          + p4*yy3(ixp,jyp,level,indexh,ngrid)
-  end do
+
+! PS clouds:
+    ip1=1
+    ip2=1
+    ip3=1
+    ip4=1
+    if (icloudbotn(ix ,jy ,indexh,ngrid) .eq. icmv) ip1=0
+    if (icloudbotn(ixp,jy ,indexh,ngrid) .eq. icmv) ip2=0
+    if (icloudbotn(ix ,jyp,indexh,ngrid) .eq. icmv) ip3=0
+    if (icloudbotn(ixp,jyp,indexh,ngrid) .eq. icmv) ip4=0
+    ipsum = ip1+ip2+ip3+ip4
+    if (ipsum .eq. 0) then
+      cbot(m)=icmv
+    else
+      cbot(m)=(ip1*p1*icloudbotn(ix ,jy ,indexh,ngrid) &
+             + ip2*p2*icloudbotn(ixp,jy ,indexh,ngrid) &
+             + ip3*p3*icloudbotn(ix ,jyp,indexh,ngrid) &
+             + ip4*p4*icloudbotn(ixp,jyp,indexh,ngrid)) / ipsum
+    endif
+
+    ip1=1
+    ip2=1
+    ip3=1
+    ip4=1
+    if (icloudthckn(ix ,jy ,indexh,ngrid) .eq. icmv) ip1=0
+    if (icloudthckn(ixp,jy ,indexh,ngrid) .eq. icmv) ip2=0
+    if (icloudthckn(ix ,jyp,indexh,ngrid) .eq. icmv) ip3=0
+    if (icloudthckn(ixp,jyp,indexh,ngrid) .eq. icmv) ip4=0
+    ipsum = ip1+ip2+ip3+ip4
+    if (ipsum .eq. 0) then
+      cthck(m)=icmv
+    else
+      cthck(m)=(ip1*p1*icloudthckn(ix ,jy ,indexh,ngrid) &
+              + ip2*p2*icloudthckn(ixp,jy ,indexh,ngrid) &
+              + ip3*p3*icloudthckn(ix ,jyp,indexh,ngrid) &
+              + ip4*p4*icloudthckn(ixp,jyp,indexh,ngrid)) / ipsum
+    endif
+! PS end clouds
+
+  enddo memind_loop
 
 
@@ -150 +198 @@
 
 
+! PS clouds:
+  icbot = nint( (cbot(1)*dt2 + cbot(2)*dt1)/dt )
+  if (nint(cbot(1)) .eq. icmv) icbot=cbot(2) 
+  if (nint(cbot(2)) .eq. icmv) icbot=cbot(1) 
+
+  icthck = nint( (cthck(1)*dt2 + cthck(2)*dt1)/dt )
+  if (nint(cthck(1)) .eq. icmv) icthck=cthck(2) 
+  if (nint(cthck(2)) .eq. icmv) icthck=cthck(1) 
+
+  if (icbot .ne. icmv .and. icthck .ne. icmv) then
+    ictop = icbot + icthck ! convert cloud thickness to cloud top
+  else
+    icbot=icmv
+    ictop=icmv
+  endif
+! PS end clouds
+
+
 end subroutine interpol_rain_nests