/[lmdze]/trunk/phylmd/readsulfate_preind.f90
ViewVC logotype

Contents of /trunk/phylmd/readsulfate_preind.f90

Parent Directory Parent Directory | Revision Log Revision Log


Revision 76 - (show annotations)
Fri Nov 15 18:45:49 2013 UTC (10 years, 6 months ago) by guez
File size: 5977 byte(s)
Moved everything out of libf.
1 SUBROUTINE readsulfate_preind (r_day, first, pi_sulfate)
2
3 ! Read in /calculate pre-industrial values of sulfate
4
5 use dimens_m
6 use dimphy
7 use temps
8 use SUPHEC_M
9 use chem
10 use getso4fromfile_m, only: getso4fromfile
11 IMPLICIT none
12
13 ! Content:
14 ! --------
15 ! This routine reads in monthly mean values of sulfate aerosols and
16 ! returns a linearly interpolated daily-mean field.
17 !
18 ! It does so for the preindustriel values of the sulfate, to a large part
19 ! analogous to the routine readsulfate.
20 !
21 ! Only Pb: Variables must be saved and don t have to be overwritten!
22 !
23 ! Author:
24 ! -------
25 ! Johannes Quaas (quaas@lmd.jussieu.fr)
26 ! 26/06/01
27 !
28 ! Input:
29 ! ------
30 double precision, intent(in):: r_day ! Day of integration
31 LOGICAL, intent(in):: first ! First timestep
32 ! (and therefore initialization necessary)
33 !
34 ! Output:
35 ! -------
36 double precision pi_sulfate (klon, klev) ! Number conc. sulfate (monthly mean data,
37 ! from file)
38 !
39 ! Local Variables:
40 ! ----------------
41 INTEGER i, ig, k, it
42 INTEGER j, iday, ny, iyr
43 parameter (ny=jjm+1)
44
45 INTEGER im, day1, day2, im2, ismaller
46 double precision pi_so4_1(iim, jjm+1, klev, 12)
47
48 double precision pi_so4(klon, klev, 12) ! SO4 in right dimension
49 SAVE pi_so4
50 double precision pi_so4_out(klon, klev)
51 SAVE pi_so4_out
52
53 CHARACTER*4 cyear
54 LOGICAL lnewday
55
56
57
58 iday = INT(r_day)
59
60 ! Get the year of the run
61 iyr = iday/360
62
63 ! Get the day of the actual year:
64 iday = iday-iyr*360
65
66 ! 0.02 is about 0.5/24, namly less than half an hour
67 lnewday = (r_day-FLOAT(iday).LT.0.02)
68
69 ! ---------------------------------------------
70 ! All has to be done only, if a new day begins!
71 ! ---------------------------------------------
72
73 IF (lnewday.OR.first) THEN
74 im = iday/30 +1 ! the actual month
75
76 ! annee_ref is the initial year (defined in temps.h)
77 iyr = iyr + annee_ref
78
79
80 IF (first) THEN
81 cyear='.nat'
82 CALL getso4fromfile(cyear,pi_so4_1)
83
84 ! Transform the horizontal 2D-field into the physics-field
85 ! (Also the levels and the latitudes have to be inversed)
86
87 ! Initialize field
88 DO it=1,12
89 DO k=1,klev
90 DO i=1,klon
91 pi_so4(i,k,it)=0.
92 ENDDO
93 ENDDO
94 ENDDO
95
96 write (*,*) 'preind: finished reading', FLOAT(iim)
97 DO it=1,12
98 DO k=1, klev
99 ! a) at the poles, use the zonal mean:
100 DO i=1,iim
101 ! North pole
102 pi_so4(1,k,it)=pi_so4(1,k,it)+pi_so4_1(i,jjm+1,klev+1-k,it)
103 ! South pole
104 pi_so4(klon,k,it)=pi_so4(klon,k,it)+pi_so4_1(i,1,klev+1-k,it)
105 ENDDO
106 pi_so4(1,k,it)=pi_so4(1,k,it)/FLOAT(iim)
107 pi_so4(klon,k,it)=pi_so4(klon,k,it)/FLOAT(iim)
108
109 ! b) the values between the poles:
110 ig=1
111 DO j=2,jjm
112 DO i=1,iim
113 ig=ig+1
114 if (ig.gt.klon) write (*,*) 'shit'
115 pi_so4(ig,k,it) = pi_so4_1(i,jjm+1-j,klev+1-k,it)
116 ENDDO
117 ENDDO
118 IF (ig.NE.klon-1) STOP 'Error in readsulfate (var conversion)'
119 ENDDO ! Loop over k (vertical)
120 ENDDO ! Loop over it (months)
121
122 ENDIF ! Had to read new data?
123
124
125 ! Interpolate to actual day:
126 IF (iday.LT.im*30-15) THEN
127 ! in the first half of the month use month before and actual month
128 im2=im-1
129 day1 = im2*30+15
130 day2 = im2*30-15
131 IF (im2.LE.0) THEN
132 ! the month is january, thus the month before december
133 im2=12
134 ENDIF
135 DO k=1,klev
136 DO i=1,klon
137 pi_sulfate(i,k) = pi_so4(i,k,im2) &
138 - FLOAT(iday-day2)/FLOAT(day1-day2) &
139 * (pi_so4(i,k,im2) - pi_so4(i,k,im))
140 IF (pi_sulfate(i,k).LT.0.) THEN
141 IF (iday-day2.LT.0.) write(*,*) 'iday-day2',iday-day2
142 IF (pi_so4(i,k,im2) - pi_so4(i,k,im).LT.0.) &
143 write(*,*) 'pi_so4(i,k,im2) - pi_so4(i,k,im)', &
144 pi_so4(i,k,im2) - pi_so4(i,k,im)
145 IF (day1-day2.LT.0.) write(*,*) 'day1-day2',day1-day2
146 stop 'pi_sulfate'
147 endif
148 ENDDO
149 ENDDO
150 ELSE
151 ! the second half of the month
152 im2=im+1
153 day1 = im*30+15
154 IF (im2.GT.12) THEN
155 ! the month is december, the following thus january
156 im2=1
157 ENDIF
158 day2 = im*30-15
159
160 DO k=1,klev
161 DO i=1,klon
162 pi_sulfate(i,k) = pi_so4(i,k,im2) &
163 - FLOAT(iday-day2)/FLOAT(day1-day2) &
164 * (pi_so4(i,k,im2) - pi_so4(i,k,im))
165 IF (pi_sulfate(i,k).LT.0.) THEN
166 IF (iday-day2.LT.0.) write(*,*) 'iday-day2',iday-day2
167 IF (pi_so4(i,k,im2) - pi_so4(i,k,im).LT.0.) &
168 write(*,*) 'pi_so4(i,k,im2) - pi_so4(i,k,im)', &
169 pi_so4(i,k,im2) - pi_so4(i,k,im)
170 IF (day1-day2.LT.0.) write(*,*) 'day1-day2',day1-day2
171 stop 'pi_sulfate'
172 endif
173 ENDDO
174 ENDDO
175 ENDIF
176
177
178 !JLD ! The sulfate concentration [molec cm-3] is read in.
179 !JLD ! Convert it into mass [ug SO4/m3]
180 !JLD ! masse_so4 in [g/mol], n_avogadro in [molec/mol]
181 DO k=1,klev
182 DO i=1,klon
183 !JLD pi_sulfate(i,k) = pi_sulfate(i,k)*masse_so4
184 !JLD . /n_avogadro*1.e+12
185 pi_so4_out(i,k) = pi_sulfate(i,k)
186 ENDDO
187 ENDDO
188
189 ELSE ! If no new day, use old data:
190 DO k=1,klev
191 DO i=1,klon
192 pi_sulfate(i,k) = pi_so4_out(i,k)
193 ENDDO
194 ENDDO
195 ENDIF ! Was this the beginning of a new day?
196
197 END SUBROUTINE readsulfate_preind

  ViewVC Help
Powered by ViewVC 1.1.21