1 |
guez |
123 |
program test_ozonecm |
2 |
|
|
|
3 |
|
|
! This is a program in Fortran 95. |
4 |
|
|
|
5 |
|
|
! This program computes values of ozone abundance from Royer, on a |
6 |
|
|
! latitude-pressure grid, and writes the values to a NetCDF file. |
7 |
|
|
! The pressure grid is "presnivs" from "disvert". |
8 |
|
|
|
9 |
|
|
use dimens_m, only: jjm, llm |
10 |
|
|
USE dimphy, ONLY : klon |
11 |
|
|
use disvert_m, only: pa, disvert, ap, bp, preff, presnivs |
12 |
|
|
use jumble, only: new_unit |
13 |
|
|
use ozonecm_m, only: ozonecm |
14 |
|
|
use phyetat0_m, only: rlat |
15 |
|
|
use nr_util, only: arth, assert |
16 |
|
|
use netcdf95, only: nf95_create, nf95_def_dim, nf95_def_var, nf95_put_att, & |
17 |
|
|
nf95_enddef, nf95_put_var, nf95_close |
18 |
|
|
use netcdf, only: nf90_clobber, nf90_float, nf90_global |
19 |
|
|
use unit_nml_m, only: unit_nml |
20 |
|
|
|
21 |
|
|
implicit none |
22 |
|
|
|
23 |
|
|
real p(llm + 1) ! pressure at LMDZ layer interface, in Pa |
24 |
|
|
real wo(klon, llm, 360) ! column density of ozone in a cell, in kDU |
25 |
|
|
real tro3(klon, llm, 360) ! mole fraction of ozone |
26 |
|
|
integer julien, k |
27 |
|
|
real, parameter:: RG = 9.80665 ! acceleration of gravity, in m s-2 |
28 |
|
|
real, parameter:: dobson_u = 2.1415e-05 ! Dobson unit, in kg m-2 |
29 |
|
|
|
30 |
|
|
! For NetCDF: |
31 |
|
|
integer ncid, dimid_time, dimid_plev, dimid_latitude |
32 |
|
|
integer varid_time, varid_plev, varid_latitude, varid_tro3 |
33 |
|
|
|
34 |
|
|
!--------------------- |
35 |
|
|
|
36 |
|
|
call assert(klon == jjm + 1, "test_ozonecm: iim should be 1") |
37 |
|
|
|
38 |
|
|
call new_unit(unit_nml) |
39 |
|
|
open(unit_nml, file="used_namelists.txt", status="replace", action="write") |
40 |
|
|
|
41 |
|
|
pa = 5e4 |
42 |
|
|
call disvert |
43 |
|
|
p = ap + bp * preff |
44 |
|
|
rlat = arth(-90., 180. / jjm, jjm + 1) |
45 |
|
|
|
46 |
|
|
do julien = 1, 360 |
47 |
|
|
wo(:, :, julien) = ozonecm(REAL(julien), spread(p, dim=1, ncopies=jjm+1)) |
48 |
|
|
end do |
49 |
|
|
|
50 |
|
|
close(unit_nml) |
51 |
|
|
|
52 |
|
|
forall (k=1: llm) |
53 |
|
|
tro3(:, k, :) & |
54 |
|
|
= wo(:, k, :) * 1e3 * dobson_u / (p(k)-p(k+1)) * rg / 48. * 29. |
55 |
|
|
end forall |
56 |
|
|
|
57 |
|
|
call nf95_create("O3_Royer.nc", nf90_clobber, ncid) |
58 |
|
|
call nf95_def_dim(ncid, "time", 360, dimid_time) |
59 |
|
|
call nf95_def_dim(ncid, "plev", llm, dimid_plev) |
60 |
|
|
call nf95_def_dim(ncid, "latitude", jjm + 1, dimid_latitude) |
61 |
|
|
call nf95_def_var(ncid, "time", nf90_float, dimid_time, varid_time) |
62 |
|
|
call nf95_def_var(ncid, "plev", nf90_float, dimid_plev, varid_plev) |
63 |
|
|
call nf95_def_var(ncid, "latitude", nf90_float, dimid_latitude, & |
64 |
|
|
varid_latitude) |
65 |
|
|
call nf95_def_var(ncid, "tro3", nf90_float, & |
66 |
|
|
(/dimid_latitude, dimid_plev, dimid_time/), varid_tro3) |
67 |
|
|
|
68 |
|
|
call nf95_put_att(ncid, varid_time, "units", "days since 1976-1-1") |
69 |
|
|
call nf95_put_att(ncid, varid_time, "calendar", "360_day") |
70 |
|
|
call nf95_put_att(ncid, varid_time, "standard_name", "time") |
71 |
|
|
|
72 |
|
|
call nf95_put_att(ncid, varid_plev, "units", "mbar") |
73 |
|
|
call nf95_put_att(ncid, varid_plev, "long_name", "pressure") |
74 |
|
|
call nf95_put_att(ncid, varid_plev, "standard_name", "air_pressure") |
75 |
|
|
|
76 |
|
|
call nf95_put_att(ncid, varid_latitude, "units", "degrees_north") |
77 |
|
|
call nf95_put_att(ncid, varid_latitude, "standard_name", "latitude") |
78 |
|
|
|
79 |
|
|
call nf95_put_att(ncid, varid_tro3, "long_name", "O3 mole fraction") |
80 |
|
|
call nf95_put_att(ncid, varid_tro3, "standard_name", & |
81 |
|
|
"mole_fraction_of_ozone_in_air") |
82 |
|
|
|
83 |
|
|
call nf95_put_att(ncid, nf90_global, "title", & |
84 |
|
|
"ozone climatology from J.-F. Royer") |
85 |
|
|
|
86 |
|
|
call nf95_enddef(ncid) |
87 |
|
|
|
88 |
|
|
call nf95_put_var(ncid, varid_time, (/(julien - 0.5, julien = 1, 360)/)) |
89 |
|
|
|
90 |
|
|
call nf95_put_var(ncid, varid_plev, presnivs(llm:1:-1) / 100.) |
91 |
|
|
! sort in ascending order and convert to hPa |
92 |
|
|
|
93 |
|
|
call nf95_put_var(ncid, varid_latitude, rlat) |
94 |
|
|
|
95 |
|
|
call nf95_put_var(ncid, varid_tro3, tro3(:, llm:1:-1, :)) |
96 |
|
|
! "plev" is in ascending order whereas "p" is descending |
97 |
|
|
|
98 |
|
|
call nf95_close(ncid) |
99 |
|
|
|
100 |
|
|
print *, 'The file "O3_Royer.nc" has been created.' |
101 |
|
|
|
102 |
|
|
end program test_ozonecm |