source: trunk/src/readdepo.f90 @ 28

Last change on this file since 28 was 4, checked in by mlanger, 11 years ago
  • Property svn:executable set to *
File size: 5.9 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 readdepo
23
24  !*****************************************************************************
25  !                                                                            *
26  !  Reads dry deposition parameters needed by the procedure of Wesely (1989). *
27  !  Wesely (1989): Parameterization of surface resistances to gaseous         *
28  !  dry deposition in regional-scale numerical models.                        *
29  !  Atmos. Environ. 23, 1293-1304.                                            *
30  !                                                                            *
31  !                                                                            *
32  !      AUTHOR: Andreas Stohl, 19 May 1995                                    *
33  !                                                                            *
34  !*****************************************************************************
35  !                                                                            *
36  ! Variables:                                                                 *
37  !                                                                            *
38  ! rcl(maxspec,5,9) [s/m] Lower canopy resistance                             *
39  ! rgs(maxspec,5,9) [s/m] Ground resistance                                   *
40  ! rlu(maxspec,5,9) [s/m] Leaf cuticular resistance                           *
41  ! rm(maxspec) [s/m]      Mesophyll resistance, set in readreleases           *
42  ! ri(maxspec) [s/m]      Stomatal resistance                                 *
43  !                                                                            *
44  ! Constants:                                                                 *
45  !                                                                            *
46  !*****************************************************************************
47
48  use par_mod
49  use com_mod
50
51  implicit none
52
53  ! FOR THIS SUBROUTINE, numclass=9 IS ASSUMED
54  !*******************************************
55
56  real :: rluh(5,numclass),rgssh(5,numclass),rgsoh(5,numclass)
57  real :: rclsh(5,numclass),rcloh(5,numclass)
58  integer :: i,j,ic
59
60
61  ! Read deposition constants related with landuse and seasonal category
62  !*********************************************************************
63  open(unitwesely,file=path(1)(1:length(1))//'surfdepo.t', &
64       status='old',err=999)
65
66  do i=1,16
67    read(unitwesely,*)
68  end do
69  do i=1,5
70    read(unitwesely,*)
71    read(unitwesely,'(8x,13f8.0)') (ri(i,j),j=1,numclass)
72    read(unitwesely,'(8x,13f8.0)') (rluh(i,j),j=1,numclass)
73    read(unitwesely,'(8x,13f8.0)') (rac(i,j),j=1,numclass)
74    read(unitwesely,'(8x,13f8.0)') (rgssh(i,j),j=1,numclass)
75    read(unitwesely,'(8x,13f8.0)') (rgsoh(i,j),j=1,numclass)
76    read(unitwesely,'(8x,13f8.0)') (rclsh(i,j),j=1,numclass)
77    read(unitwesely,'(8x,13f8.0)') (rcloh(i,j),j=1,numclass)
78  end do
79
80  ! TEST
81  ! do 31 i=1,5
82  !   ri(i,13)=ri(i,5)
83  !   rluh(i,13)=rluh(i,5)
84  !   rac(i,13)=rac(i,5)
85  !   rgssh(i,13)=rgssh(i,5)
86  !   rgsoh(i,13)=rgsoh(i,5)
87  !   rclsh(i,13)=rclsh(i,5)
88  !   rcloh(i,13)=rcloh(i,5)
89  !31             continue
90  ! TEST
91  ! Sabine Eckhardt, Dec 06, set resistances of 9999 to 'infinite' (1E25)
92  do i=1,5
93    do j=1,numclass
94      if    (ri(i,j).eq.9999.)    ri(i,j)=1.E25
95      if  (rluh(i,j).eq.9999.)  rluh(i,j)=1.E25
96      if   (rac(i,j).eq.9999.)   rac(i,j)=1.E25
97      if (rgssh(i,j).eq.9999.) rgssh(i,j)=1.E25
98      if (rgsoh(i,j).eq.9999.) rgsoh(i,j)=1.E25
99      if (rclsh(i,j).eq.9999.) rclsh(i,j)=1.E25
100      if (rcloh(i,j).eq.9999.) rcloh(i,j)=1.E25
101    end do
102  end do
103
104
105
106  do i=1,5
107    do j=1,numclass
108      ri(i,j)=max(ri(i,j),0.001)
109      rluh(i,j)=max(rluh(i,j),0.001)
110      rac(i,j)=max(rac(i,j),0.001)
111      rgssh(i,j)=max(rgssh(i,j),0.001)
112      rgsoh(i,j)=max(rgsoh(i,j),0.001)
113      rclsh(i,j)=max(rclsh(i,j),0.001)
114      rcloh(i,j)=max(rcloh(i,j),0.001)
115    end do
116  end do
117  close(unitwesely)
118
119
120  ! Compute additional parameters
121  !******************************
122
123  do ic=1,nspec
124    if (reldiff(ic).gt.0.) then     ! gas is dry deposited
125      do i=1,5
126        do j=1,numclass
127          rlu(ic,i,j)=rluh(i,j)/(1.e-5*henry(ic)+f0(ic))
128          rgs(ic,i,j)=1./(henry(ic)/(10.e5*rgssh(i,j))+f0(ic)/ &
129               rgsoh(i,j))
130          rcl(ic,i,j)=1./(henry(ic)/(10.e5*rclsh(i,j))+f0(ic)/ &
131               rcloh(i,j))
132        end do
133      end do
134    endif
135  end do
136
137
138  return
139
140
141999   write(*,*) '### FLEXPART ERROR! FILE              ###'
142  write(*,*) '### surfdepo.t DOES NOT EXIST.        ###'
143  stop
144
145end subroutine readdepo
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG