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