/[MITgcm]/MITgcm/pkg/fizhi/update_chemistry_exports.F
ViewVC logotype

Contents of /MITgcm/pkg/fizhi/update_chemistry_exports.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.11 - (show annotations) (download)
Fri Oct 22 14:52:14 2004 UTC (19 years, 8 months ago) by molod
Branch: MAIN
CVS Tags: checkpoint57d_post, checkpoint57b_post, checkpoint57c_pre, checkpoint55j_post, checkpoint56b_post, checkpoint57e_post, checkpoint56c_post, checkpoint57a_post, checkpoint57a_pre, checkpoint55i_post, checkpoint57, checkpoint56, eckpoint57e_pre, checkpoint57c_post, checkpoint56a_post
Changes since 1.10: +3 -2 lines
Change myTime from integer to real (_RL), change pressure0 calculation in ini_fixed
(Little bugs.....)

1 C $Header: /u/gcmpack/MITgcm/pkg/fizhi/update_chemistry_exports.F,v 1.10 2004/07/28 01:25:07 molod Exp $
2 C $Name: $
3
4 #include "FIZHI_OPTIONS.h"
5 subroutine update_chemistry_exports (myTime, myIter, myThid)
6 c----------------------------------------------------------------------
7 c Subroutine update_chemistry_exports - 'Wrapper' routine to update
8 c the fields related to the earth's chemistry that are needed
9 c by fizhi.
10 c Also: Set up "bi, bj loop" and some timers and clocks here.
11 c
12 c Call: interp_chemistry
13 c-----------------------------------------------------------------------
14 implicit none
15 #include "SIZE.h"
16 #include "fizhi_SIZE.h"
17 #include "fizhi_land_SIZE.h"
18 #include "GRID.h"
19 #include "DYNVARS.h"
20 #include "fizhi_chemistry_coms.h"
21 #include "fizhi_coms.h"
22 #include "gridalt_mapping.h"
23 #include "EEPARAMS.h"
24 #include "chronos.h"
25
26 integer myIter, myThid
27 _RL myTime
28
29 c pe on physics grid refers to bottom edge
30 _RL pephy(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nrphys+1,nSx,nSy)
31 _RL pphy(sNx,sNy,Nrphys,nSx,nSy)
32 _RL oz1(nlatsoz,nlevsoz), strq1(nlatsq,nlevsq)
33 _RL waterin(sNx,sNy,Nrphys), xlat(sNx,sNy)
34
35 integer i, j, L, LL, bi, bj
36 integer im1, im2, jm1, jm2
37 integer nhms1,nymd1,nhms2,nymd2,imns,ipls
38 _RL facm, facp
39 logical alarm
40 external alarm
41
42 im1 = 1
43 im2 = sNx
44 jm1 = 1
45 jm2 = sNy
46
47 if( alarm('radsw').or.alarm('radlw') ) then
48
49 do bj = myByLo(myThid), myByHi(myThid)
50 do bi = myBxLo(myThid), myBxHi(myThid)
51
52 c Construct the physics grid pressures - count pephy levels top down
53 c (even though dpphy counted bottom up)
54 do j = 1,sNy
55 do i = 1,sNx
56 pephy(i,j,Nrphys+1,bi,bj)=(Ro_surf(i,j,bi,bj)+etaH(i,j,bi,bj))
57 do L = 2,Nrphys+1
58 LL = Nrphys+2-L
59 pephy(i,j,LL,bi,bj)=pephy(i,j,LL+1,bi,bj)-dpphys(i,j,L-1,bi,bj)
60 enddo
61 enddo
62 enddo
63 do j = 1,sNy
64 do i = 1,sNx
65 do L = 1,Nrphys
66 pphy(i,j,L,bi,bj)=(pephy(i,j,L+1,bi,bj)+pephy(i,j,L,bi,bj))
67 . /200.
68 enddo
69 enddo
70 enddo
71
72 do j = 1,sNy
73 do i = 1,sNx
74 xlat(i,j) = yC(i,j,bi,bj)
75 do L = 1,Nrphys
76 waterin(i,j,L) = sphy(i,j,L,bi,bj)
77 enddo
78 enddo
79 enddo
80
81 call time_bound(nymd,nhms,nymd1,nhms1,nymd2,nhms2,imns,ipls)
82 call interp_time(nymd,nhms,nymd1,nhms1,nymd2,nhms2,facm,facp)
83
84 do L = 1,nlevsoz
85 do j = 1,nlatsoz
86 oz1(j,L) = ozone(j,L,imns)*facm + ozone(j,L,ipls)*facp
87 enddo
88 enddo
89
90 do L = 1,nlevsq
91 do j = 1,nlatsq
92 strq1(j,L) = stratq(j,L,imns)*facm + stratq(j,L,ipls)*facp
93 enddo
94 enddo
95
96 call interp_chemistry(strq1,nlevsq,nlatsq,levsq,latsq,
97 . oz1,nlevsoz,nlatsoz,levsoz,latsoz,waterin,pphy,xlat,
98 . im2,jm2,Nrphys,nSx,nSy,bi,bj,o3,qstr)
99
100 enddo
101 enddo
102
103 endif
104
105 return
106 end
107
108 subroutine interp_chemistry (stratq,nwatlevs,nwatlats,watlevs,
109 . watlats,ozone,nozlevs,nozlats,ozlevs,ozlats,
110 . qz,plz,xlat,im,jm,lm,nSx,nSy,bi,bj,ozrad,qzrad)
111
112 implicit none
113
114 c Input Variables
115 c ---------------
116 integer nwatlevs,nwatlats,nozlevs,nozlats,nSx,nSy,bi,bj
117 _RL stratq(nwatlats,nwatlevs),ozone(nozlats,nozlevs)
118 _RL watlevs(nwatlevs),watlats(nwatlats)
119 _RL ozlevs(nozlevs),ozlats(nozlats)
120 integer im,jm,lm
121 _RL qz(im,jm,lm),plz(im,jm,lm)
122 _RL xlat(im,jm)
123 _RL ozrad(im,jm,lm,nSx,nSy)
124 _RL qzrad(im,jm,lm,nSx,nSy)
125
126 C **********************************************************************
127 C **** Get Ozone and Stratospheric Moisture Data ****
128 C **********************************************************************
129
130 call interp_qz (stratq,nwatlevs,nwatlats,watlevs,watlats,im*jm,
131 . bi,bj, xlat,lm,plz,qz,qzrad(1,1,1,bi,bj))
132 call interp_oz (ozone ,nozlevs,nozlats,ozlevs,ozlats,im*jm,
133 . bi,bj, xlat,lm,plz,ozrad(1,1,1,bi,bj))
134
135 return
136 end
137
138 subroutine interp_qz(stratq,nwatlevs,nwatlats,watlevs,watlats,
139 . irun,bi,bj,xlat,nlevs,pres,qz_in,qz_out )
140 C***********************************************************************
141 C Purpose
142 C To Interpolate Chemistry Moisture from Chemistry Grid to Physics Grid
143 C
144 C INPUT Argument Description
145 C stratq .... Climatological SAGE Stratospheric Moisture
146 C irun ...... Number of Columns to be filled
147 C xlat ...... Latitude in Degrees
148 C nlevs ..... Vertical Dimension
149 C pres ...... PRES (IM,JM,nlevs) Three-dimensional array of pressures
150 C qz_in ..... Model Moisture (kg/kg mass mixing radtio)
151 C qz_out .... Combination of Chemistry Moisture and Model Moisture
152 C (kg/kg mass mixing ratio)
153 C
154 C***********************************************************************
155
156 implicit none
157 integer nwatlevs,nwatlats
158 integer bi,bj
159 _RL stratq ( nwatlats,nwatlevs )
160 _RL watlats (nwatlats)
161 _RL watlevs (nwatlevs)
162
163 integer irun,nlevs
164 _RL xlat (irun)
165 _RL pres (irun,nlevs)
166 _RL qz_in (irun,nlevs)
167 _RL qz_out(irun,nlevs)
168
169 c Local Variables
170 c ---------------
171 integer pqu,pql,dpq
172 parameter ( pqu = 100. )
173 parameter ( pql = 300. )
174 parameter ( dpq = pql-pqu )
175
176 integer i,k,L1,L2,LM,LP
177 _RL h2o_time_lat (irun,nwatlevs)
178 _RL qz_clim(irun,nlevs)
179
180 _RL qpr1(irun), qpr2(irun), slope(irun)
181 _RL pr1(irun), pr2(irun)
182
183 integer jlat,jlatm,jlatp
184
185 C **********************************************************************
186 C **** Interpolate Moisture data to model latitudes ***
187 C **********************************************************************
188
189 DO 32 k = 1, nwatlevs
190 DO 34 i = 1,irun
191
192 DO 36 jlat = 1, nwatlats
193 IF( watlats(jlat).gt.xlat(i) ) THEN
194 IF( jlat.EQ.1 ) THEN
195 jlatm = 1
196 jlatp = 1
197 slope(i) = 0
198 ELSE
199 jlatm = jlat -1
200 jlatp = jlat
201 slope(i) = ( xlat(i) -watlats(jlat-1) )
202 . / ( watlats(jlat)-watlats(jlat-1) )
203 ENDIF
204 GOTO 37
205 ENDIF
206 36 CONTINUE
207 jlatm = nwatlats
208 jlatp = nwatlats
209 slope(i) = 1
210 37 CONTINUE
211 QPR1(i) = stratq(jlatm,k)
212 QPR2(i) = stratq(jlatp,k)
213 34 CONTINUE
214
215 do i = 1,irun
216 h2o_time_lat(i,k) = qpr1(i) + slope(i)*(qpr2(i)-qpr1(i))
217 enddo
218
219 32 CONTINUE
220
221 C **********************************************************************
222 C **** Interpolate Latitude Moisture data to model pressures ***
223 C **********************************************************************
224
225 DO 40 L2 = 1,nlevs
226
227 DO 44 i= 1, irun
228 DO 46 L1 = 1,nwatlevs
229 IF( watlevs(L1).GT.pres(i,L2) ) THEN
230 IF( L1.EQ.1 ) THEN
231 LM = 1
232 LP = 2
233 ELSE
234 LM = L1-1
235 LP = L1
236 ENDIF
237 GOTO 47
238 ENDIF
239 46 CONTINUE
240 LM = nwatlevs-1
241 LP = nwatlevs
242 47 CONTINUE
243 PR1(i) = watlevs (LM)
244 PR2(i) = watlevs (LP)
245 QPR1(i) = h2o_time_lat(i,LM)
246 QPR2(i) = h2o_time_lat(i,LP)
247 44 CONTINUE
248
249 do i= 1, irun
250 slope(i) =(QPR1(i)-QPR2(i)) / (PR1(i)-PR2(i))
251 qz_clim(i,L2) = QPR2(i) + (pres(i,L2)-PR2(i))*SLOPE(i)
252 enddo
253
254 40 CONTINUE
255
256 c
257 c ... Above 100 mb, using climatological water data set ...................
258 c ... Below 300 mb, using model predicted water data set ...................
259 c ... In between, using linear interpolation ...............................
260 c
261 do k= 1, nlevs
262 do i= 1, irun
263 if( pres(i,k).ge.pqu .and. pres(i,k).le. pql) then
264 qz_out(i,k) = qz_clim(i,k)+(qz_in(i,k)-
265 1 qz_clim(i,k))*(pres(i,k)-pqu)/dpq
266 else if( pres(i,k) .gt. pql ) then
267 qz_out(i,k) = qz_in (i,k)
268 else
269 qz_out(i,k) = qz_clim(i,k)
270 endif
271 enddo
272 enddo
273
274 return
275 end
276
277 subroutine interp_oz (ozone,nozlevs,nozlats,ozlevs,ozlats,
278 . irun,bi,bj,xlat,nlevs,plevs,ozrad)
279 C***********************************************************************
280 C Purpose
281 C To Interpolate Chemistry Ozone from Chemistry Grid to Physics Grid
282 C
283 C INPUT Argument Description
284 C ozone ..... Climatological Ozone
285 C chemistry .. Chemistry State Data Structure
286 C irun ....... Number of Columns to be filled
287 C xlat ....... Latitude in Degrees
288 C nlevs ...... Vertical Dimension
289 C pres ....... Three-dimensional array of pressures
290 C ozrad ...... Ozone on Physics Grid (kg/kg mass mixing radtio)
291 C
292 C***********************************************************************
293 implicit none
294 integer nozlevs,nozlats,irun,nlevs
295 integer bi,bj
296 _RL ozone(nozlats,nozlevs)
297 _RL xlat(irun)
298 _RL plevs(irun,nlevs)
299 _RL ozrad(irun,nlevs)
300 _RL ozlevs(nozlevs),ozlats(nozlats)
301
302 c Local Variables
303 c ---------------
304 _RL zero,one,o3min,voltomas
305 PARAMETER ( ZERO = 0.0 )
306 PARAMETER ( ONE = 1.0 )
307 PARAMETER ( O3MIN = 1.0E-10 )
308 PARAMETER ( VOLTOMAS = 1.655E-6 )
309
310 integer i,k,L1,L2,LM,LP
311 integer jlat,jlatm,jlatp
312 _RL O3INT1(IRUN,nozlevs)
313 _RL QPR1(IRUN), QPR2(IRUN), SLOPE(IRUN)
314 _RL PR1(IRUN), PR2(IRUN)
315
316 C **********************************************************************
317 C **** INTERPOLATE ozone data to model latitudes ***
318 C **********************************************************************
319
320 DO 32 K=1,nozlevs
321 DO 34 I=1,IRUN
322
323 DO 36 jlat = 1,nozlats
324 IF( ozlats(jlat).gt.xlat(i) ) THEN
325 IF( jlat.EQ.1 ) THEN
326 jlatm = 1
327 jlatp = 1
328 slope(i) = zero
329 ELSE
330 jlatm = jlat-1
331 jlatp = jlat
332 slope(i) = ( XLAT(I) -ozlats(jlat-1) )
333 . / ( ozlats(jlat)-ozlats(jlat-1) )
334 ENDIF
335 GOTO 37
336 ENDIF
337 36 CONTINUE
338 jlatm = nozlats
339 jlatp = nozlats
340 slope(i) = one
341 37 CONTINUE
342 QPR1(I) = ozone(jlatm,k)
343 QPR2(I) = ozone(jlatp,k)
344 34 CONTINUE
345
346 DO 38 I=1,IRUN
347 o3int1(i,k) = qpr1(i) + slope(i)*( qpr2(i)-qpr1(i) )
348 38 CONTINUE
349
350 32 CONTINUE
351
352 C **********************************************************************
353 C **** INTERPOLATE latitude ozone data to model pressures ***
354 C **********************************************************************
355
356 DO 40 L2 = 1,NLEVS
357
358 DO 44 I = 1,IRUN
359 DO 46 L1 = 1,nozlevs
360 IF( ozlevs(L1).GT.PLEVS(I,L2) ) THEN
361 IF( L1.EQ.1 ) THEN
362 LM = 1
363 LP = 2
364 ELSE
365 LM = L1-1
366 LP = L1
367 ENDIF
368 GOTO 47
369 ENDIF
370 46 CONTINUE
371 LM = nozlevs-1
372 LP = nozlevs
373 47 CONTINUE
374 PR1(I) = ozlevs (LM)
375 PR2(I) = ozlevs (LP)
376 QPR1(I) = O3INT1(I,LM)
377 QPR2(I) = O3INT1(I,LP)
378 44 CONTINUE
379
380 DO 48 I=1,IRUN
381 SLOPE(I) = ( QPR1(I)-QPR2(I) )
382 . / ( PR1(I)- PR2(I) )
383 ozrad(I,L2) = QPR2(I) + ( PLEVS(I,L2)-PR2(I) )*SLOPE(I)
384
385 if( ozrad(i,l2).lt.o3min ) then
386 ozrad(i,l2) = o3min
387 endif
388
389 48 CONTINUE
390 40 CONTINUE
391
392 C **********************************************************************
393 C **** CONVERT FROM VOLUME MIXING RATIO TO MASS MIXING RATIO ***
394 C **********************************************************************
395
396 DO 60 I=1,IRUN*NLEVS
397 ozrad (I,1) = ozrad(I,1) * VOLTOMAS
398 60 CONTINUE
399
400 RETURN
401 END
402

  ViewVC Help
Powered by ViewVC 1.1.22