/[lmdze]/trunk/libf/phylmd/ozonecm.f
ViewVC logotype

Contents of /trunk/libf/phylmd/ozonecm.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3 - (show annotations)
Wed Feb 27 13:16:39 2008 UTC (16 years, 2 months ago) by guez
File size: 3224 byte(s)
Initial import
1 !
2 ! $Header: /home/cvsroot/LMDZ4/libf/phylmd/ozonecm.F,v 1.3 2005/06/06 13:16:33 fairhead Exp $
3 !
4 SUBROUTINE ozonecm(rjour, rlat, paprs, o3)
5 use dimens_m
6 use dimphy
7 use clesphys
8 use YOMCST
9 IMPLICIT none
10 C
11 C The ozone climatology is based on an analytic formula which fits the
12 C Krueger and Mintzner (1976) profile, as well as the variations with
13 C altitude and latitude of the maximum ozone concentrations and the total
14 C column ozone concentration of Keating and Young (1986). The analytic
15 C formula have been established by J-F Royer (CRNM, Meteo France), who
16 C also provided us the code.
17 C
18 C A. J. Krueger and R. A. Minzner, A Mid-Latitude Ozone Model for the
19 C 1976 U.S. Standard Atmosphere, J. Geophys. Res., 81, 4477, (1976).
20 C
21 C Keating, G. M. and D. F. Young, 1985: Interim reference models for the
22 C middle atmosphere, Handbook for MAP, vol. 16, 205-229.
23 C
24
25 real, intent(in):: rjour
26 REAL, intent(in):: rlat(klon)
27 real, intent(in):: paprs(klon,klev+1)
28 REAL, intent(out):: o3(klon,klev) ! ozone concentration in kg/kg
29
30 REAL tozon
31 real pi, pl
32 INTEGER i, k
33 C----------------------------------------------------------
34 REAL field(klon,klev+1)
35 REAL ps
36 PARAMETER (ps=101325.0)
37 REAL an, unit, zo3q3
38 SAVE an, unit, zo3q3
39 REAL mu,gms, zslat, zsint, zcost, z, ppm, qpm, a
40 REAL asec, bsec, aprim, zo3a3
41 C----------------------------------------------------------
42 c data an /365.25/ (meteo)
43 DATA an /360.00/
44 DATA unit /2.1415e-05/
45 DATA zo3q3 /4.0e-08/
46
47 pi = 4.0 * ATAN(1.0)
48 DO k = 1, klev
49 DO i = 1, klon
50 zslat = SIN(pi/180.*rlat(i))
51 zsint = SIN(2.*pi*(rjour+15.)/an)
52 zcost = COS(2.*pi*(rjour+15.)/an)
53 z = 0.0531+zsint*(-0.001595+0.009443*zslat) +
54 . zcost*(-0.001344-0.00346*zslat) +
55 . zslat**2*(.056222+zslat**2*(-.037609
56 . +.012248*zsint+.00521*zcost+.008890*zslat))
57 zo3a3 = zo3q3/ps/2.
58 z = z-zo3q3*ps
59 gms = z
60 mu = ABS(sin(pi/180.*rlat(i)))
61 ppm = 800.-(500.*zslat+150.*zcost)*zslat
62 qpm = 1.74e-5-(7.5e-6*zslat+1.7e-6*zcost)*zslat
63 bsec = 2650.+5000.*zslat**2
64 a = 4.0*(bsec)**(3./2.)*(ppm)**(3./2.)*(1.0+(bsec/ps)**(3./2.))
65 a = a/(bsec**(3./2.)+ppm**(3./2.))**2
66 aprim = (2.666666*qpm*ppm-a*gms)/(1.0-a)
67 aprim = amax1(0.,aprim)
68 asec = (gms-aprim)*(1.0+(bsec/ps)**(3./2.))
69 asec = AMAX1(0.0,asec)
70 aprim = gms-asec/(1.+(bsec/ps)**(3./2.))
71 pl = paprs(i,k)
72 tozon = aprim/(1.+3.*(ppm/pl)**2)+asec/(1.+(bsec/pl)**(3./2.))
73 . + zo3a3*pl*pl
74 tozon = tozon / 9.81 ! en kg/m**2
75 tozon = tozon / unit ! en kg/m**2 > u dobson (10e-5 m)
76 tozon = tozon / 1000. ! en cm
77 field(i,k) = tozon
78 ENDDO
79 ENDDO
80 c
81 DO i = 1, klon
82 field(i,klev+1) = 0.0
83 ENDDO
84 DO k = 1, klev
85 DO i = 1, klon
86 o3(i,k) = field(i,k) - field(i,k+1)
87 IF (.not. bug_ozone) then
88 c convert o3 into kg/kg
89 o3(i,k)=MAX(o3(i,k),1.0e-12)*RG/46.6968
90 . /(paprs(i,k)-paprs(i,k+1))
91 ENDIF
92 ENDDO
93 ENDDO
94 c
95 RETURN
96 END

  ViewVC Help
Powered by ViewVC 1.1.21