phylmd/Mobidic/regr_pr_int.f

module regr_pr_int_m

  ! Author: Lionel GUEZ

  implicit none

contains

  subroutine regr_pr_int(ncid, name, julien, pplay, top_value, v3)

    ! "regr_pr_int" stands for "regrid pressure interpolation".

    ! This procedure reads a 2D latitude-pressure field from a NetCDF
    ! file, at a given day, regrids the field in pressure to the LMDZ
    ! vertical grid and packs it to the LMDZ horizontal "physics"
    ! grid. We assume that, in the input file, the field has 3
    ! dimensions: latitude, pressure, julian day. We assume that
    ! latitudes are in ascending order in the input file. We assume
    ! that the input data is already on the LMDZ "rlatu" latitude
    ! grid.

    ! The target vertical LMDZ grid is the grid of mid-layers. The
    ! input data does not depend on longitude, but the pressure at
    ! LMDZ mid-layers does. Therefore, the values on the LMDZ grid do
    ! depend on longitude. Regridding is by linear interpolation.

    use dimens_m, only: iim, jjm, llm
    use dimphy, only: klon
    use grid_change, only: gr_dyn_phy
    use netcdf95, only: nf95_inq_varid, nf95_get_var
    use nr_util, only: assert
    use numer_rec_95, only: regr1_lint
    use press_coefoz_m, only: plev

    integer, intent(in):: ncid ! NetCDF ID of the file
    character(len=*), intent(in):: name ! of the NetCDF variable
    integer, intent(in):: julien ! jour julien, 1 <= julien <= 360

    real, intent(in):: pplay(:, :) ! (klon, llm)
    ! (pression pour le mileu de chaque couche, en Pa)

    real, intent(in):: top_value
    ! (extra value of field at 0 pressure)

    real, intent(out):: v3(:, :) ! (klon, llm)
    ! Regridded field on the partial "physics" grid. "v3(i, k)" is at
    ! longitude "xlon(i)", latitude "xlat(i)", middle of layer "k".

    ! Variables local to the procedure:

    integer varid ! for NetCDF

    real v1(jjm + 1, 0:size(plev))
    ! Input field at day "julien". "v1(j, k >=1)" is at latitude
    ! "rlatu(j)" and pressure "plev(k)".

    real v2(klon, 0:size(plev))
    ! Field on the "physics" horizontal grid.  "v2(i, k >= 1)" is at
    ! longitude "xlon(i)", latitude "xlat(i)" and pressure "plev(k)".)

    integer i

    !--------------------------------------------

    call assert(shape(v3) == (/klon, llm/), "regr_pr_int")

    call nf95_inq_varid(ncid, name, varid)

    ! Get data at the right day from the input file:
    call nf95_get_var(ncid, varid, v1(:, 1:), start=(/1, 1, julien/))
    ! Latitudes are in ascending order in the input file while
    ! "rlatu" is in descending order so we need to invert order:
    v1(:, 1:) = v1(jjm+1:1:-1, 1:)

    ! Complete "v1" with the value at 0 pressure:
    v1(:, 0) = top_value

    v2 = gr_dyn_phy(spread(v1, dim = 1, ncopies = iim + 1))

    ! Regrid in pressure at each horizontal position:
    do i = 1, klon
       v3(i, llm:1:-1) = regr1_lint(v2(i, :), (/0., plev/), pplay(i, llm:1:-1))
       ! (invert order of indices because "pplay" is in descending order)
    end do

  end subroutine regr_pr_int

end module regr_pr_int_m
1	guez	90	module regr_pr_int_m
2
3	guez	168	! Author: Lionel GUEZ
4
5	guez	90	implicit none
6
7			contains
8
9	guez	162	subroutine regr_pr_int(ncid, name, julien, pplay, top_value, v3)
10	guez	90
11			! "regr_pr_int" stands for "regrid pressure interpolation".
12
13			! This procedure reads a 2D latitude-pressure field from a NetCDF
14			! file, at a given day, regrids the field in pressure to the LMDZ
15			! vertical grid and packs it to the LMDZ horizontal "physics"
16			! grid. We assume that, in the input file, the field has 3
17			! dimensions: latitude, pressure, julian day. We assume that
18			! latitudes are in ascending order in the input file. We assume
19			! that the input data is already on the LMDZ "rlatu" latitude
20			! grid.
21
22	guez	173	! The target vertical LMDZ grid is the grid of mid-layers. The
23			! input data does not depend on longitude, but the pressure at
24			! LMDZ mid-layers does. Therefore, the values on the LMDZ grid do
25			! depend on longitude. Regridding is by linear interpolation.
26	guez	90
27			use dimens_m, only: iim, jjm, llm
28			use dimphy, only: klon
29	guez	168	use grid_change, only: gr_dyn_phy
30	guez	90	use netcdf95, only: nf95_inq_varid, nf95_get_var
31			use nr_util, only: assert
32			use numer_rec_95, only: regr1_lint
33			use press_coefoz_m, only: plev
34
35			integer, intent(in):: ncid ! NetCDF ID of the file
36			character(len=*), intent(in):: name ! of the NetCDF variable
37			integer, intent(in):: julien ! jour julien, 1 <= julien <= 360
38
39	guez	162	real, intent(in):: pplay(:, :) ! (klon, llm)
40			! (pression pour le mileu de chaque couche, en Pa)
41
42	guez	90	real, intent(in):: top_value
43			! (extra value of field at 0 pressure)
44
45			real, intent(out):: v3(:, :) ! (klon, llm)
46			! Regridded field on the partial "physics" grid. "v3(i, k)" is at
47			! longitude "xlon(i)", latitude "xlat(i)", middle of layer "k".
48
49			! Variables local to the procedure:
50
51	guez	105	integer varid ! for NetCDF
52	guez	90
53	guez	162	real v1(jjm + 1, 0:size(plev))
54	guez	90	! Input field at day "julien". "v1(j, k >=1)" is at latitude
55			! "rlatu(j)" and pressure "plev(k)".
56
57	guez	162	real v2(klon, 0:size(plev))
58			! Field on the "physics" horizontal grid. "v2(i, k >= 1)" is at
59			! longitude "xlon(i)", latitude "xlat(i)" and pressure "plev(k)".)
60	guez	90
61	guez	168	integer i
62	guez	90
63			!--------------------------------------------
64
65			call assert(shape(v3) == (/klon, llm/), "regr_pr_int")
66
67			call nf95_inq_varid(ncid, name, varid)
68
69			! Get data at the right day from the input file:
70			call nf95_get_var(ncid, varid, v1(:, 1:), start=(/1, 1, julien/))
71			! Latitudes are in ascending order in the input file while
72			! "rlatu" is in descending order so we need to invert order:
73			v1(:, 1:) = v1(jjm+1:1:-1, 1:)
74
75			! Complete "v1" with the value at 0 pressure:
76			v1(:, 0) = top_value
77
78	guez	168	v2 = gr_dyn_phy(spread(v1, dim = 1, ncopies = iim + 1))
79	guez	162
80	guez	90	! Regrid in pressure at each horizontal position:
81	guez	162	do i = 1, klon
82			v3(i, llm:1:-1) = regr1_lint(v2(i, :), (/0., plev/), pplay(i, llm:1:-1))
83			! (invert order of indices because "pplay" is in descending order)
84	guez	90	end do
85
86			end subroutine regr_pr_int
87
88			end module regr_pr_int_m