Changeset 37 for branches/petra/src/interpol_rain.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.f90 (modified) (7 diffs)

branches/petra/src/interpol_rain.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   *
@@ -20 +20 @@
 !**********************************************************************
 
-subroutine interpol_rain(yy1,yy2,yy3,nxmax,nymax,nzmax,nx, &
-       ny,memind,xt,yt,level,itime1,itime2,itime,yint1,yint2,yint3)
+subroutine interpol_rain(yy1,yy2,yy3,nxmax,nymax,nzmax,nx,ny, &
+       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     o     o     o
@@ -40 +41 @@
   !                                                                           *
   !     30 August 1996                                                        *
+  !****************************************************************************
+  ! Petra Seibert, 2011/2012-2015:
+  !  Add interpolation of cloud bottom and cloud thickness
+  !  for fix to SE's new wet scavenging scheme
   !                                                                           *
   !****************************************************************************
@@ -66 +71 @@
   !                                                                           *
   !****************************************************************************
+  
+  use com_mod, only: icloudbot, icloudthck
+  use par_mod, only: icmv
 
   implicit none
 
   integer :: nx,ny,nxmax,nymax,nzmax,memind(2),m,ix,jy,ixp,jyp
-  integer :: itime,itime1,itime2,level,indexh
+  integer :: itime,itime1,itime2,level,indexh,ip1,ip2,ip3,ip4 
+  integer :: icbot,icthck,ictop,ipsum
   real :: yy1(0:nxmax-1,0:nymax-1,nzmax,2)
   real :: yy2(0:nxmax-1,0:nymax-1,nzmax,2)
   real :: yy3(0:nxmax-1,0:nymax-1,nzmax,2)
-  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
-
 
 
@@ -112 +120 @@
   !***********************
 
+  memind_loop: &
   do m=1,2
+
     indexh=memind(m)
 
@@ -127 +137 @@
          + p3*yy3(ix ,jyp,level,indexh) &
          + p4*yy3(ixp,jyp,level,indexh)
-  end do
+
+
+! PS clouds:
+    ip1=1
+    ip2=1
+    ip3=1
+    ip4=1
+    if (icloudbot(ix ,jy ,indexh) .eq. icmv) ip1=0
+    if (icloudbot(ixp,jy ,indexh) .eq. icmv) ip2=0
+    if (icloudbot(ix ,jyp,indexh) .eq. icmv) ip3=0
+    if (icloudbot(ixp,jyp,indexh) .eq. icmv) ip4=0
+    ipsum = ip1+ip2+ip3+ip4
+    if (ipsum .eq. 0) then
+      cbot(m)=icmv
+    else
+      cbot(m)=(ip1*p1*icloudbot(ix ,jy ,indexh) &
+             + ip2*p2*icloudbot(ixp,jy ,indexh) &
+             + ip3*p3*icloudbot(ix ,jyp,indexh) &
+             + ip4*p4*icloudbot(ixp,jyp,indexh)) / ipsum
+    endif
+
+    ip1=1
+    ip2=1
+    ip3=1
+    ip4=1
+    if (icloudthck(ix ,jy ,indexh) .eq. icmv) ip1=0
+    if (icloudthck(ixp,jy ,indexh) .eq. icmv) ip2=0
+    if (icloudthck(ix ,jyp,indexh) .eq. icmv) ip3=0
+    if (icloudthck(ixp,jyp,indexh) .eq. icmv) ip4=0
+    ipsum = ip1+ip2+ip3+ip4
+    if (ipsum .eq. 0) then
+      cthck(m)=icmv
+    else
+      cthck(m)=(ip1*p1*icloudthck(ix ,jy ,indexh) &
+              + ip2*p2*icloudthck(ixp,jy ,indexh) &
+              + ip3*p3*icloudthck(ix ,jyp,indexh) &
+              + ip4*p4*icloudthck(ixp,jyp,indexh)) / ipsum
+    endif
+! PS end clouds
+
+  enddo memind_loop
 
 
@@ -142 +192 @@
   yint3=(y3(1)*dt2+y3(2)*dt1)/dt
 
+! 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