source: flexpart.git/src/interpol_rain_nests.f90 @ 02095e3

10.4.1_peseiGFS_025bugfixes+enhancementsdevrelease-10release-10.4.1scaling-bugunivie
Last change on this file since 02095e3 was db712a8, checked in by Espen Sollum ATMOS <eso@…>, 8 years ago

Completed handling of nested wind fields with cloud water (for wet deposition).

  • Property mode set to 100644
File size: 7.6 KB
Line 
1!**********************************************************************
2! Copyright 1998,1999,2000,2001,2002,2005,2007,2008,2009,2010         *
3! Andreas Stohl, Petra Seibert, A. Frank, Gerhard Wotawa,             *
4! Caroline Forster, Sabine Eckhardt, John Burkhart, Harald Sodemann   *
5!                                                                     *
6! This file is part of FLEXPART.                                      *
7!                                                                     *
8! FLEXPART is free software: you can redistribute it and/or modify    *
9! it under the terms of the GNU General Public License as published by*
10! the Free Software Foundation, either version 3 of the License, or   *
11! (at your option) any later version.                                 *
12!                                                                     *
13! FLEXPART is distributed in the hope that it will be useful,         *
14! but WITHOUT ANY WARRANTY; without even the implied warranty of      *
15! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the       *
16! GNU General Public License for more details.                        *
17!                                                                     *
18! You should have received a copy of the GNU General Public License   *
19! along with FLEXPART.  If not, see <http://www.gnu.org/licenses/>.   *
20!**********************************************************************
21
22subroutine interpol_rain_nests(yy1,yy2,yy3,nxmaxn,nymaxn,nzmax, &
23       maxnests,ngrid,nxn,nyn,memind,xt,yt,level,itime1,itime2,itime, &
24       yint1,yint2,yint3)
25  !                                i   i   i    i      i      i
26  !   i       i    i   i    i    i  i    i     i      i      i
27  !  o     o     o
28  !****************************************************************************
29  !                                                                           *
30  !  Interpolation of meteorological fields on 2-d model layers for nested    *
31  !  grids. This routine is related to levlin3interpol.f for the mother domain*
32  !                                                                           *
33  !  In horizontal direction bilinear interpolation interpolation is used.    *
34  !  Temporally a linear interpolation is used.                               *
35  !  Three fields are interpolated at the same time.                          *
36  !                                                                           *
37  !  This is a special version of levlininterpol to save CPU time.            *
38  !                                                                           *
39  !  1 first time                                                             *
40  !  2 second time                                                            *
41  !                                                                           *
42  !                                                                           *
43  !     Author: A. Stohl                                                      *
44  !                                                                           *
45  !     15 March 2000                                                         *
46  !                                                                           *
47  !****************************************************************************
48  !                                                                           *
49  ! Variables:                                                                *
50  !                                                                           *
51  ! dt1,dt2              time differences between fields and current position *
52  ! dz1,dz2              z distance between levels and current position       *
53  ! height(nzmax)        heights of the model levels                          *
54  ! indexh               help variable                                        *
55  ! indz                 the level closest to the current trajectory position *
56  ! indzh                help variable                                        *
57  ! itime                current time                                         *
58  ! itime1               time of the first wind field                         *
59  ! itime2               time of the second wind field                        *
60  ! ix,jy                x,y coordinates of lower left subgrid point          *
61  ! level                level at which interpolation shall be done           *
62  ! memind(3)            points to the places of the wind fields              *
63  ! nx,ny                actual field dimensions in x,y and z direction       *
64  ! nxmax,nymax,nzmax    maximum field dimensions in x,y and z direction      *
65  ! xt                   current x coordinate                                 *
66  ! yint                 the final interpolated value                         *
67  ! yt                   current y coordinate                                 *
68  ! yy(0:nxmax,0:nymax,nzmax,3) meteorological field used for interpolation   *
69  ! zt                   current z coordinate                                 *
70  !                                                                           *
71  !****************************************************************************
72  use par_mod, only: numwfmem
73
74  implicit none
75
76  integer :: maxnests,ngrid
77  integer :: nxn(maxnests),nyn(maxnests),nxmaxn,nymaxn,nzmax,memind(numwfmem)
78  integer :: m,ix,jy,ixp,jyp,itime,itime1,itime2,level,indexh
79  real :: yy1(0:nxmaxn-1,0:nymaxn-1,nzmax,numwfmem,maxnests)
80  real :: yy2(0:nxmaxn-1,0:nymaxn-1,nzmax,numwfmem,maxnests)
81  real :: yy3(0:nxmaxn-1,0:nymaxn-1,nzmax,numwfmem,maxnests)
82  real :: ddx,ddy,rddx,rddy,dt1,dt2,dt,y1(2),y2(2),y3(2)
83  real :: xt,yt,yint1,yint2,yint3,p1,p2,p3,p4
84
85
86
87  ! If point at border of grid -> small displacement into grid
88  !***********************************************************
89
90  ! if (xt.ge.(real(nxn(ngrid)-1)-0.0001)) &
91  !      xt=real(nxn(ngrid)-1)-0.0001
92  ! if (yt.ge.(real(nyn(ngrid)-1)-0.0001)) &
93  !      yt=real(nyn(ngrid)-1)-0.0001
94
95! ESO make it consistent with interpol_rain
96  if (xt.ge.(real(nxn(ngrid)-1))) xt=real(nxn(ngrid)-1)-0.00001
97  if (yt.ge.(real(nyn(ngrid)-1))) yt=real(nyn(ngrid)-1)-0.00001
98
99
100
101  !**********************************************************************
102  ! 1.) Bilinear horizontal interpolation
103  ! This has to be done separately for 2 fields (Temporal)
104  !*******************************************************
105
106  ! Determine the lower left corner and its distance to the current position
107  !*************************************************************************
108
109  ix=int(xt)
110  jy=int(yt)
111
112  ixp=ix+1
113  jyp=jy+1
114  ddx=xt-real(ix)
115  ddy=yt-real(jy)
116  rddx=1.-ddx
117  rddy=1.-ddy
118  p1=rddx*rddy
119  p2=ddx*rddy
120  p3=rddx*ddy
121  p4=ddx*ddy
122
123
124  ! Loop over 2 time steps
125  !***********************
126
127  do m=1,2
128    indexh=memind(m)
129
130    y1(m)=p1*yy1(ix ,jy ,level,indexh,ngrid) &
131         + p2*yy1(ixp,jy ,level,indexh,ngrid) &
132         + p3*yy1(ix ,jyp,level,indexh,ngrid) &
133         + p4*yy1(ixp,jyp,level,indexh,ngrid)
134    y2(m)=p1*yy2(ix ,jy ,level,indexh,ngrid) &
135         + p2*yy2(ixp,jy ,level,indexh,ngrid) &
136         + p3*yy2(ix ,jyp,level,indexh,ngrid) &
137         + p4*yy2(ixp,jyp,level,indexh,ngrid)
138    y3(m)=p1*yy3(ix ,jy ,level,indexh,ngrid) &
139         + p2*yy3(ixp,jy ,level,indexh,ngrid) &
140         + p3*yy3(ix ,jyp,level,indexh,ngrid) &
141         + p4*yy3(ixp,jyp,level,indexh,ngrid)
142  end do
143
144
145  !************************************
146  ! 2.) Temporal interpolation (linear)
147  !************************************
148
149  dt1=real(itime-itime1)
150  dt2=real(itime2-itime)
151  dt=dt1+dt2
152
153  yint1=(y1(1)*dt2+y1(2)*dt1)/dt
154  yint2=(y2(1)*dt2+y2(2)*dt1)/dt
155  yint3=(y3(1)*dt2+y3(2)*dt1)/dt
156
157
158end subroutine interpol_rain_nests
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG