1 | ! |
---|
2 | ! $Id: o3cm.F 1299 2010-01-20 14:27:21Z jghattas $ |
---|
3 | ! |
---|
4 | SUBROUTINE o3cm (amb, bmb, sortie, ntab) |
---|
5 | IMPLICIT none |
---|
6 | c====================================================================== |
---|
7 | c Auteur(s): Z.X. Li (LMD/CNRS) date: 19930818 |
---|
8 | c Objet: Ce programme calcule le contenu en ozone "sortie" |
---|
9 | c (unite: cm.atm) entre deux niveaux "amb" et "bmb" (unite: mb) |
---|
10 | c "ntab" est le nombre d'intervalles pour l'integration, sa |
---|
11 | c valeur depend bien sur de l'epaisseur de la couche et de |
---|
12 | c la precision qu'on souhaite a obtenir |
---|
13 | c====================================================================== |
---|
14 | REAL amb, bmb, sortie |
---|
15 | INTEGER ntab |
---|
16 | c====================================================================== |
---|
17 | INTEGER n |
---|
18 | REAL xtab(500), xa, xb, ya, yb, xincr |
---|
19 | c====================================================================== |
---|
20 | external mbtozm |
---|
21 | CHARACTER (LEN=20) :: modname='' |
---|
22 | CHARACTER (LEN=80) :: abort_message |
---|
23 | c====================================================================== |
---|
24 | c la fonction en ligne w(x) donne le profil de l'ozone en fonction |
---|
25 | c de l'altitude (unite: cm.atm / km) |
---|
26 | c (Green 1964, Appl. Opt. 3: 203-208) |
---|
27 | REAL wp, xp, h, x, w, con |
---|
28 | PARAMETER (wp=0.218, xp=23.25, h=4.63, con=1.0) |
---|
29 | w(x) = wp/h * EXP((x-xp)/h)/ (con+EXP((x-xp)/h))**2 |
---|
30 | c====================================================================== |
---|
31 | IF (ntab .GT. 499) THEN |
---|
32 | abort_message = 'BIG ntab' |
---|
33 | CALL abort_gcm (modname,abort_message,1) |
---|
34 | ENDIF |
---|
35 | xincr = (bmb-amb) / REAL(ntab) |
---|
36 | xtab(1) = amb |
---|
37 | DO n = 2, ntab |
---|
38 | xtab(n) = xtab(n-1) + xincr |
---|
39 | ENDDO |
---|
40 | xtab(ntab+1) = bmb |
---|
41 | sortie = 0.0 |
---|
42 | DO n = 1, ntab |
---|
43 | CALL mbtozm(xtab(n), xa) |
---|
44 | CALL mbtozm(xtab(n+1), xb) |
---|
45 | xa = xa / 1000. |
---|
46 | xb = xb / 1000. |
---|
47 | ya = w(xa) |
---|
48 | yb = w(xb) |
---|
49 | sortie = sortie + (ya+yb)/2.0 * ABS(xb-xa) |
---|
50 | ENDDO |
---|
51 | RETURN |
---|
52 | END |
---|
53 | SUBROUTINE mbtozm(rmb,zm) |
---|
54 | IMPLICIT none |
---|
55 | c====================================================================== |
---|
56 | c Auteur(s): Z.X. Li (LMD/CNRS) |
---|
57 | c Objet: transformer une hauteur de mb (rmb) en metre (zm) |
---|
58 | c====================================================================== |
---|
59 | REAL rmb, zm |
---|
60 | c====================================================================== |
---|
61 | REAL gama, tzero, pzero, g, r |
---|
62 | PARAMETER (gama=6.5e-3, tzero=288., pzero=1013.25) |
---|
63 | PARAMETER (g=9.81, r=287.0) |
---|
64 | zm = tzero/gama * ( 1.-(rmb/pzero)**(r*gama/g) ) |
---|
65 | RETURN |
---|
66 | END |
---|