[524] | 1 | ! |
---|
[1299] | 2 | ! $Id: o3cm.F 1299 2010-01-20 14:27:21Z fhourdin $ |
---|
[524] | 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 |
---|
[1299] | 21 | CHARACTER (LEN=20) :: modname='' |
---|
| 22 | CHARACTER (LEN=80) :: abort_message |
---|
[524] | 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====================================================================== |
---|
[1299] | 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) |
---|
[524] | 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 |
---|