/[lmdze]/trunk/Sources/phylmd/readsulfate.f
ViewVC logotype

Contents of /trunk/Sources/phylmd/readsulfate.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 191 - (show annotations)
Mon May 9 19:56:28 2016 UTC (8 years ago) by guez
File size: 7753 byte(s)
Extracted the call to read_comdissnew out of conf_gcm.

Made ok_instan a variable of module clesphys, itau_phy a variable of
module phyetat0_m, nid_ins a variable of module ini_histins_m, itap a
variable of new module time_phylmdz, so that histwrite_phy can be
called from any procedure without the need to cascade those variables
into that procedure. Made itau_w a variable of module time_phylmdz so
that it is computed only once per time step of physics.

Extracted variables of module clesphys which were in namelist
conf_phys_nml into their own namelist, clesphys_nml, and created
procedure read_clesphys reading clesphys_nml, to avoid side effect.

No need for double precision in procedure getso4fromfile. Assume there
is a single variable for the whole year in the NetCDF file instead of
one variable per month.

Created generic procedure histwrite_phy and removed procedure
write_histins, following LMDZ. histwrite_phy has only two arguments,
can be called from anywhere, and should manage the logic of writing or
not writing into various history files with various operations. So the
test on ok_instan goes inside histwrite_phy.

Test for raz_date in phyetat0 instead of physiq to avoid side effect.

Created procedure increment_itap to avoid side effect.

Removed unnecessary differences between procedures readsulfate and
readsulfate_pi.

1 module readsulfate_m
2
3 IMPLICIT none
4
5 contains
6
7 SUBROUTINE readsulfate(dayvrai, time, first, sulfate)
8
9 ! From LMDZ4/libf/phylmd/readsulfate.F, version 1.2 2005/05/19
10 ! 08:27:15 fairhead
11
12 ! This routine reads monthly mean values of sulfate aerosols and
13 ! returns a linearly interpolated daily-mean field.
14
15 ! Author: Johannes Quaas (quaas@lmd.jussieu.fr)
16 ! April 26th, 2001
17
18 ! ATTENTION!! runs are supposed to start with Jan, 1. 1930
19 ! (rday = 1)
20
21 ! The model always has 360 days per year.
22 ! SO4 concentration rather than mixing ratio
23 ! 10yr-mean-values to interpolate
24 ! Introduce flag to read in just one decade
25
26 use dimens_m, only: iim, jjm
27 use dimphy, only: klon, klev
28 use dynetat0_m, only: annee_ref
29 use getso4fromfile_m, only: getso4fromfile
30
31 integer, intent(in):: dayvrai
32 ! current day number, based at value 1 on January 1st of annee_ref
33
34 REAL, intent(in):: time ! heure de la journ\'ee en fraction de jour
35
36 LOGICAL, intent(in):: first ! First timestep
37 ! (and therefore initialization necessary)
38
39 real, intent(out):: sulfate(klon, klev)
40 ! mass of sulfate (monthly mean data, from file) (micro g SO4 / m3)
41
42 ! Local:
43 INTEGER i, ig, k, it
44 INTEGER j, iday, iyr, iyr1, iyr2
45 CHARACTER(len = 4) cyear
46 INTEGER im, day1, day2, im2
47 real so4_1(iim, jjm + 1, klev, 12)
48 real so4_2(iim, jjm + 1, klev, 12) ! sulfate distributions
49 double precision, save:: so4(klon, klev, 12) ! SO4 in right dimension
50 double precision, save:: so4_out(klon, klev)
51 LOGICAL lnewday
52
53 !----------------------------------------------------------------
54
55 iday = dayvrai
56
57 ! Get the year of the run
58 iyr = iday/360
59
60 ! Get the day of the actual year:
61 iday = iday - iyr*360
62
63 ! 0.02 is about 0.5/24, namly less than half an hour
64 lnewday = time < 0.02
65
66 ! All has to be done only, if a new day begins
67
68 test_newday: IF (lnewday .OR. first) THEN
69 im = iday/30 + 1 ! the actual month
70
71 ! annee_ref is the initial year (defined in temps.h)
72 iyr = iyr + annee_ref
73
74 ! Do I have to read new data? (Is this the first day of a year?)
75 IF (first .OR. iday == 1.) THEN
76 ! Initialize field
77 DO it = 1, 12
78 DO k = 1, klev
79 DO i = 1, klon
80 so4(i, k, it) = 0.
81 ENDDO
82 ENDDO
83 ENDDO
84
85 IF (iyr < 1850) THEN
86 cyear = '.nat'
87 print *, 'getso4 iyr = ', iyr, ' ', cyear
88 CALL getso4fromfile(cyear, so4_1)
89 ELSE IF (iyr >= 2100) THEN
90 cyear = '2100'
91 print *, 'getso4 iyr = ', iyr, ' ', cyear
92 CALL getso4fromfile(cyear, so4_1)
93 ELSE
94 ! Read in data:
95 ! a) from actual 10-yr-period
96
97 IF (iyr < 1900) THEN
98 iyr1 = 1850
99 iyr2 = 1900
100 ELSE IF (iyr >= 1900.and.iyr < 1920) THEN
101 iyr1 = 1900
102 iyr2 = 1920
103 ELSE
104 iyr1 = INT(iyr/10)*10
105 iyr2 = INT(1 + iyr/10)*10
106 ENDIF
107 WRITE(cyear, '(I4)') iyr1
108 print *, 'getso4 iyr = ', iyr, ' ', cyear
109 CALL getso4fromfile(cyear, so4_1)
110
111 ! Read two decades:
112 ! b) from the next following one
113 WRITE(cyear, '(I4)') iyr2
114 print *, 'getso4 iyr = ', iyr, ' ', cyear
115 CALL getso4fromfile(cyear, so4_2)
116
117 ! Interpolate linarily to the actual year:
118 DO it = 1, 12
119 DO k = 1, klev
120 DO j = 1, jjm
121 DO i = 1, iim
122 so4_1(i, j, k, it) = so4_1(i, j, k, it) &
123 - REAL(iyr - iyr1)/REAL(iyr2 - iyr1) &
124 * (so4_1(i, j, k, it) - so4_2(i, j, k, it))
125 ENDDO
126 ENDDO
127 ENDDO
128 ENDDO
129 ENDIF
130
131 ! Transform the horizontal 2D-field into the physics-field
132 ! (Also the levels and the latitudes have to be inversed)
133
134 DO it = 1, 12
135 DO k = 1, klev
136 ! a) at the poles, use the zonal mean:
137 DO i = 1, iim
138 ! North pole
139 so4(1, k, it) = so4(1, k, it) &
140 + so4_1(i, jjm + 1, klev + 1 - k, it)
141 ! South pole
142 so4(klon, k, it) = so4(klon, k, it) &
143 + so4_1(i, 1, klev + 1 - k, it)
144 ENDDO
145 so4(1, k, it) = so4(1, k, it)/REAL(iim)
146 so4(klon, k, it) = so4(klon, k, it)/REAL(iim)
147
148 ! b) the values between the poles:
149 ig = 1
150 DO j = 2, jjm
151 DO i = 1, iim
152 ig = ig + 1
153 if (ig > klon) stop 1
154 so4(ig, k, it) = so4_1(i, jjm + 1 - j, klev + 1 - k, it)
155 ENDDO
156 ENDDO
157 IF (ig /= klon - 1) then
158 print *, 'Error in readsulfate (var conversion)'
159 STOP 1
160 end IF
161 ENDDO ! Loop over k (vertical)
162 ENDDO ! Loop over it (months)
163 ENDIF ! Had to read new data?
164
165 ! Interpolate to actual day:
166 IF (iday < im*30 - 15) THEN
167 ! in the first half of the month use month before and actual month
168 im2 = im - 1
169 day1 = im2*30 + 15
170 day2 = im2*30 - 15
171 IF (im2 <= 0) THEN
172 ! the month is january, thus the month before december
173 im2 = 12
174 ENDIF
175 DO k = 1, klev
176 DO i = 1, klon
177 sulfate(i, k) = so4(i, k, im2) &
178 - REAL(iday - day2)/REAL(day1 - day2) &
179 * (so4(i, k, im2) - so4(i, k, im))
180 IF (sulfate(i, k) < 0.) THEN
181 IF (iday - day2 < 0.) write(*, *) 'iday - day2', iday - day2
182 IF (so4(i, k, im2) - so4(i, k, im) < 0.) &
183 write(*, *) 'so4(i, k, im2) - so4(i, k, im)', &
184 so4(i, k, im2) - so4(i, k, im)
185 IF (day1 - day2 < 0.) write(*, *) 'day1 - day2', day1 - day2
186 stop 1
187 endif
188 ENDDO
189 ENDDO
190 ELSE
191 ! the second half of the month
192 im2 = im + 1
193 IF (im2 > 12) THEN
194 ! the month is december, the following thus january
195 im2 = 1
196 ENDIF
197 day2 = im*30 - 15
198 day1 = im*30 + 15
199 DO k = 1, klev
200 DO i = 1, klon
201 sulfate(i, k) = so4(i, k, im2) &
202 - REAL(iday - day2)/REAL(day1 - day2) &
203 * (so4(i, k, im2) - so4(i, k, im))
204 IF (sulfate(i, k) < 0.) THEN
205 IF (iday - day2 < 0.) write(*, *) 'iday - day2', iday - day2
206 IF (so4(i, k, im2) - so4(i, k, im) < 0.) &
207 write(*, *) 'so4(i, k, im2) - so4(i, k, im)', &
208 so4(i, k, im2) - so4(i, k, im)
209 IF (day1 - day2 < 0.) write(*, *) 'day1 - day2', day1 - day2
210 stop 1
211 endif
212 ENDDO
213 ENDDO
214 ENDIF
215
216 DO k = 1, klev
217 DO i = 1, klon
218 so4_out(i, k) = sulfate(i, k)
219 ENDDO
220 ENDDO
221 ELSE
222 ! If no new day, use old data:
223 DO k = 1, klev
224 DO i = 1, klon
225 sulfate(i, k) = so4_out(i, k)
226 ENDDO
227 ENDDO
228 ENDIF test_newday
229
230 END SUBROUTINE readsulfate
231
232 end module readsulfate_m

  ViewVC Help
Powered by ViewVC 1.1.21