source: LMDZ4/trunk/libf/phylmd/ozonecm.F @ 648

Last change on this file since 648 was 641, checked in by Laurent Fairhead, 20 years ago

Rajout d'un commentaire a la routine pour s'y retrouver MPL/JLD
LF

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 2.7 KB
Line 
1!
2! $Header$
3!
4      SUBROUTINE ozonecm(rjour, rlat, paprs, o3)
5      IMPLICIT none
6C
7C The ozone climatology is based on an analytic formula which fits the
8C Krueger and Mintzner (1976) profile, as well as the variations with
9C altitude and latitude of the maximum ozone concentrations and the total
10C column ozone concentration of Keating and Young (1986). The analytic
11C formula have been established by J-F Royer (CRNM, Meteo France), who
12C also provided us the code.
13C
14C A. J. Krueger and R. A. Minzner, A Mid-Latitude Ozone Model for the
15C 1976 U.S. Standard Atmosphere, J. Geophys. Res., 81, 4477, (1976).
16C
17C Keating, G. M. and D. F. Young, 1985: Interim reference models for the
18C middle atmosphere, Handbook for MAP, vol. 16, 205-229.
19C
20
21#include "dimensions.h"
22#include "dimphy.h"
23      REAL rlat(klon), paprs(klon,klev+1)
24      REAL o3(klon,klev)
25      REAL tozon, rjour, pi, pl
26      INTEGER i, k
27C----------------------------------------------------------
28      REAL field(klon,klev+1)
29      REAL ps
30      PARAMETER (ps=101325.0)
31      REAL an, unit, zo3q3
32      SAVE an, unit, zo3q3
33      REAL mu,gms, zslat, zsint, zcost, z, ppm, qpm, a
34      REAL asec, bsec, aprim, zo3a3
35C----------------------------------------------------------
36c         data an /365.25/   (meteo)
37      DATA an /360.00/
38      DATA unit /2.1415e-05/
39      DATA zo3q3 /4.0e-08/
40
41      pi = 4.0 * ATAN(1.0)
42      DO k = 1, klev
43      DO i = 1, klon
44      zslat = SIN(pi/180.*rlat(i))
45      zsint = SIN(2.*pi*(rjour+15.)/an)
46      zcost = COS(2.*pi*(rjour+15.)/an)
47      z = 0.0531+zsint*(-0.001595+0.009443*zslat) +
48     .  zcost*(-0.001344-0.00346*zslat) +
49     .  zslat**2*(.056222+zslat**2*(-.037609
50     . +.012248*zsint+.00521*zcost+.008890*zslat))
51      zo3a3 = zo3q3/ps/2.
52      z = z-zo3q3*ps
53      gms = z
54      mu = ABS(sin(pi/180.*rlat(i)))
55      ppm = 800.-(500.*zslat+150.*zcost)*zslat
56      qpm = 1.74e-5-(7.5e-6*zslat+1.7e-6*zcost)*zslat
57      bsec = 2650.+5000.*zslat**2
58      a = 4.0*(bsec)**(3./2.)*(ppm)**(3./2.)*(1.0+(bsec/ps)**(3./2.))
59      a = a/(bsec**(3./2.)+ppm**(3./2.))**2
60      aprim = (2.666666*qpm*ppm-a*gms)/(1.0-a)
61      aprim = amax1(0.,aprim)
62      asec = (gms-aprim)*(1.0+(bsec/ps)**(3./2.))
63      asec = AMAX1(0.0,asec)
64      aprim = gms-asec/(1.+(bsec/ps)**(3./2.))
65      pl = paprs(i,k)
66      tozon = aprim/(1.+3.*(ppm/pl)**2)+asec/(1.+(bsec/pl)**(3./2.))
67     .  + zo3a3*pl*pl
68      tozon = tozon / 9.81  ! en kg/m**2
69      tozon = tozon / unit  ! en kg/m**2  > u dobson (10e-5 m)
70      tozon = tozon / 1000. ! en cm
71      field(i,k) = tozon
72      ENDDO
73      ENDDO
74c
75      DO i = 1, klon
76         field(i,klev+1) = 0.0
77      ENDDO
78      DO k = 1, klev
79      DO i = 1, klon
80         o3(i,k) = field(i,k) - field(i,k+1)
81      ENDDO
82      ENDDO
83c
84      RETURN
85      END
Note: See TracBrowser for help on using the repository browser.