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

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

Parent Directory Parent Directory | Revision Log Revision Log

Revision 130 - (hide annotations)
Tue Feb 24 15:43:51 2015 UTC (9 years, 3 months ago) by guez
Original Path: trunk/phylmd/readsulfate.f
File size: 8438 byte(s) 
The information in argument rdayvrai of calfis was redundant with the
information in argument time. Furthermore, in the physics part of gcm,
we need separately the day number (an integer) and the time of
day. So, replaced real argument rdayvrai of calfis containing elapsed
time by integer argument dayvrai containing day number. Corresponding
change in leapfrog. In procedure physiq, replaced real argument
rdayvrai by integer argument dayvrai. In procedures readsulfate and
readsulfate_preind, replaced real argument r_day by arguments dayvrai
and time.

In procedure alboc, replaced real argument rjour by integer argument
jour. alboc was always called by interfsurf_hq with actual argument
real(jour), and the meaning of the dummy argument in alboc seems to be
that it should be an integer.

In procedure leapfrog, local variable time could not be > 1. Removed
test.

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