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

Last change on this file since 955 was 766, checked in by Laurent Fairhead, 17 years ago

Merge entre la version V3_conv et le HEAD
YM, JG, LF

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