/[lmdze]/trunk/phylmd/Mobidic/regr_pr_o3.f
ViewVC logotype

Annotation of /trunk/phylmd/Mobidic/regr_pr_o3.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 313 - (hide annotations)
Mon Dec 10 15:54:30 2018 UTC (5 years, 5 months ago) by guez
File size: 4441 byte(s)
Remove module temps. Move variable itau_dyn from module temps to
module dynetat0_m, where it is defined.

Split module dynetat0_m into dynetat0_m and dynetat0_chosen_m. The
motivation is to create smaller modules. Procedures principal_cshift
and invert_zoomx had to stay in dynetat0_m because of circular
dependency. Now we will be able to move them away. Module variables
which are chosen by the user, not computed, in program ce0l go to
dynetat0_chosen_m: day_ref, annee_ref, clon, clat, grossismx,
grossismy, dzoomx, dzoomy, taux, tauy.

Move variable "pa" from module disvert_m to module
dynetat0_chosen_m. Define "pa" in dynetat0_chosen rather than etat0.

Define day_ref and annee_ref in procedure read_serre rather than
etat0.

1 guez 8 module regr_pr_o3_m
2    
3     implicit none
4    
5     contains
6    
7 guez 97 subroutine regr_pr_o3(p3d, o3_mob_regr)
8 guez 8
9     ! "regr_pr_o3" stands for "regrid pressure ozone".
10     ! This procedure reads Mobidic ozone mole fraction from
11     ! "coefoz_LMDZ.nc" at the initial day and regrids it in pressure.
12 guez 19 ! Ozone mole fraction from "coefoz_LMDZ.nc" is a 2D latitude --
13     ! pressure variable.
14     ! The target horizontal LMDZ grid is the "scalar" grid: "rlonv", "rlatu".
15     ! The target vertical LMDZ grid is the grid of layer boundaries.
16     ! We assume that the input variable is already on the LMDZ "rlatu"
17     ! latitude grid.
18     ! The input variable does not depend on longitude, but the
19     ! pressure at LMDZ layers does.
20     ! Therefore, the values on the LMDZ grid do depend on longitude.
21 guez 8 ! Regridding is by averaging, assuming a step function.
22 guez 19 ! We assume that, in the input file, the pressure levels are in
23     ! hPa and strictly increasing.
24 guez 8
25 guez 265 use dimensions, only: iim, jjm, llm
26 guez 313 use dynetat0_chosen_m, only: day_ref
27 guez 129 use grid_change, only: dyn_phy
28 guez 313 use netcdf, only: nf90_nowrite
29     use netcdf95, only: nf95_open, nf95_close, nf95_inq_varid, nf95_get_var, &
30 guez 22 nf95_gw_var
31 guez 36 use nr_util, only: assert
32 guez 61 use numer_rec_95, only: regr1_step_av
33 guez 8
34 guez 97 REAL, intent(in):: p3d(:, :, :) ! (iim + 1, jjm + 1, llm+1)
35     ! pressure at layer interfaces, in Pa
36     ! ("p3d(i, j, l)" is at longitude "rlonv(i)", latitude "rlatu(j)",
37     ! for interface "l")
38    
39 guez 8 real, intent(out):: o3_mob_regr(:, :, :) ! (iim + 1, jjm + 1, llm)
40     ! (ozone mole fraction from Mobidic adapted to the LMDZ grid)
41     ! ("o3_mob_regr(i, j, l)" is at longitude "rlonv(i)", latitude
42     ! "rlatu(j)" and pressure level "pls(i, j, l)")
43    
44     ! Variables local to the procedure:
45    
46 guez 225 real, allocatable:: plev(:)
47 guez 8 ! (pressure levels of Mobidic data, in Pa, in strictly increasing order)
48    
49     real, allocatable:: press_in_edg(:)
50     ! (edges of pressure intervals for Mobidic data, in Pa, in strictly
51     ! increasing order)
52    
53 guez 313 integer ncid, varid ! for NetCDF
54 guez 8 integer n_plev ! number of pressure levels in Mobidic data
55 guez 19 integer i, j
56 guez 8
57     real, allocatable:: r_mob(:, :)! (jjm + 1, n_plev)
58 guez 129 ! (ozone mole fraction from Mobidic at day "day_ref")
59 guez 8 ! (r_mob(j, k) is at latitude "rlatu(j)" and pressure level "plev(k)".)
60    
61     !------------------------------------------------------------
62    
63 guez 10 print *, "Call sequence information: regr_pr_o3"
64 guez 8
65 guez 97 call assert(shape(o3_mob_regr) == (/iim + 1, jjm + 1, llm/), &
66     "regr_pr_o3 o3_mob_regr")
67     call assert(shape(p3d) == (/iim + 1, jjm + 1, llm + 1/), &
68     "regr_pr_o3 p3d")
69    
70 guez 8 call nf95_open("coefoz_LMDZ.nc", nf90_nowrite, ncid)
71    
72 guez 22 call nf95_inq_varid(ncid, "plev", varid)
73     call nf95_gw_var(ncid, varid, plev)
74 guez 8 ! Convert from hPa to Pa because "regr_pr_av" requires so:
75     plev = plev * 100.
76     n_plev = size(plev)
77    
78     ! Compute edges of pressure intervals:
79     allocate(press_in_edg(n_plev + 1))
80     press_in_edg(1) = 0.
81     ! We choose edges halfway in logarithm:
82     forall (j = 2:n_plev) press_in_edg(j) = sqrt(plev(j - 1) * plev(j))
83     press_in_edg(n_plev + 1) = huge(0.)
84     ! (infinity, but any value guaranteed to be greater than the
85     ! surface pressure would do)
86    
87     call nf95_inq_varid(ncid, "r_Mob", varid)
88     allocate(r_mob(jjm + 1, n_plev))
89    
90     ! Get data at the right day from the input file:
91 guez 313 call nf95_get_var(ncid, varid, r_mob, start=(/1, 1, day_ref/))
92 guez 8 ! Latitudes are in increasing order in the input file while
93     ! "rlatu" is in decreasing order so we need to invert order:
94     r_mob = r_mob(jjm+1:1:-1, :)
95    
96     call nf95_close(ncid)
97    
98 guez 20 ! Regrid in pressure by averaging a step function of pressure:
99 guez 19 do j = 1, jjm + 1
100     do i = 1, iim
101     if (dyn_phy(i, j)) then
102     o3_mob_regr(i, j, llm:1:-1) &
103     = regr1_step_av(r_mob(j, :), press_in_edg, &
104     p3d(i, j, llm+1:1:-1))
105     ! (invert order of indices because "p3d" is decreasing)
106     end if
107     end do
108     end do
109 guez 8
110 guez 19 ! Duplicate pole values on all longitudes:
111     o3_mob_regr(2:, 1, :) = spread(o3_mob_regr(1, 1, :), dim=1, ncopies=iim)
112     o3_mob_regr(2:, jjm + 1, :) &
113     = spread(o3_mob_regr(1, jjm + 1, :), dim=1, ncopies=iim)
114    
115     ! Duplicate first longitude to last longitude:
116     o3_mob_regr(iim + 1, 2:jjm, :) = o3_mob_regr(1, 2:jjm, :)
117    
118 guez 8 end subroutine regr_pr_o3
119    
120     end module regr_pr_o3_m

  ViewVC Help
Powered by ViewVC 1.1.21