[1158] | 1 | ! $Id$ |
---|
| 2 | module regr_pr |
---|
| 3 | |
---|
| 4 | ! In both procedures of this module: |
---|
| 5 | ! -- the root process reads a 2D latitude-pressure field from a |
---|
| 6 | ! NetCDF file, at a given day. |
---|
| 7 | ! -- the field is packed to the LMDZ horizontal "physics" |
---|
| 8 | ! grid and scattered to all threads of all processes; |
---|
| 9 | ! -- in all the threads of all the processes, the field is regridded in |
---|
| 10 | ! pressure to the LMDZ vertical grid. |
---|
| 11 | ! We assume that, in the input file, the field has 3 dimensions: |
---|
| 12 | ! latitude, pressure, julian day. |
---|
| 13 | ! We assume that latitudes are in ascending order in the input file. |
---|
| 14 | |
---|
| 15 | |
---|
| 16 | implicit none |
---|
| 17 | |
---|
| 18 | contains |
---|
| 19 | |
---|
| 20 | subroutine regr_pr_av(ncid, name, julien, press_in_edg, paprs, v3) |
---|
| 21 | |
---|
| 22 | ! "regr_pr_av" stands for "regrid pressure averaging". |
---|
| 23 | ! The target vertical LMDZ grid is the grid of layer boundaries. |
---|
| 24 | ! Regridding in pressure is done by averaging a step function of pressure. |
---|
| 25 | |
---|
| 26 | use dimphy, only: klon |
---|
| 27 | use netcdf95, only: nf95_inq_varid, handle_err |
---|
| 28 | use netcdf, only: nf90_get_var |
---|
| 29 | use assert_m, only: assert |
---|
| 30 | use regr1_step_av_m, only: regr1_step_av |
---|
| 31 | use mod_phys_lmdz_mpi_data, only: is_mpi_root |
---|
| 32 | |
---|
| 33 | use mod_phys_lmdz_transfert_para, only: scatter2d |
---|
| 34 | ! (pack to the LMDZ horizontal "physics" grid and scatter) |
---|
| 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 | real, intent(in):: press_in_edg(:) |
---|
| 41 | ! (edges of pressure intervals for input data, in Pa, in strictly |
---|
| 42 | ! ascending order) |
---|
| 43 | |
---|
| 44 | real, intent(in):: paprs(:, :) ! (klon, llm + 1) |
---|
| 45 | ! (pression pour chaque inter-couche, en Pa) |
---|
| 46 | |
---|
| 47 | real, intent(out):: v3(:, :) ! (klon, llm) |
---|
| 48 | ! (regridded field on the partial "physics" grid) |
---|
| 49 | ! ("v3(i, k)" is at longitude "xlon(i)", latitude |
---|
| 50 | ! "xlat(i)", in pressure interval "[paprs(i, l+1), paprs(i, l)]".) |
---|
| 51 | |
---|
| 52 | ! Variables local to the procedure: |
---|
| 53 | |
---|
| 54 | include "dimensions.h" |
---|
| 55 | integer varid, ncerr ! for NetCDF |
---|
| 56 | |
---|
| 57 | real v1(iim, jjm + 1, size(press_in_edg) - 1) |
---|
| 58 | ! (input field at day "julien", on the global "dynamics" horizontal grid) |
---|
| 59 | ! (First dimension is for longitude. |
---|
| 60 | ! The value is the same for all longitudes. |
---|
| 61 | ! "v1(:, j, k)" is at latitude "rlatu(j)" and for |
---|
| 62 | ! pressure interval "[press_in_edg(k), press_in_edg(k+1)]".) |
---|
| 63 | |
---|
| 64 | real v2(klon, size(press_in_edg) - 1) |
---|
| 65 | ! (field scattered to the partial "physics" horizontal grid) |
---|
| 66 | ! ("v2(i, k)" is at longitude "xlon(i)", latitude "xlat(i)" and |
---|
| 67 | ! for pressure interval "[press_in_edg(k), press_in_edg(k+1)]".) |
---|
| 68 | |
---|
| 69 | integer i |
---|
| 70 | |
---|
| 71 | !-------------------------------------------- |
---|
| 72 | |
---|
| 73 | call assert(shape(v3) == (/klon, llm/), "regr_pr_av v3") |
---|
| 74 | call assert(shape(paprs) == (/klon, llm+1/), "regr_pr_av paprs") |
---|
| 75 | |
---|
| 76 | !$omp master |
---|
| 77 | if (is_mpi_root) then |
---|
| 78 | call nf95_inq_varid(ncid, name, varid) |
---|
| 79 | |
---|
| 80 | ! Get data at the right day from the input file: |
---|
| 81 | ncerr = nf90_get_var(ncid, varid, v1(1, :, :), start=(/1, 1, julien/)) |
---|
| 82 | call handle_err("regr_pr_av nf90_get_var " // name, ncerr, ncid) |
---|
| 83 | ! Latitudes are in ascending order in the input file while |
---|
| 84 | ! "rlatu" is in descending order so we need to invert order: |
---|
| 85 | v1(1, :, :) = v1(1, jjm+1:1:-1, :) |
---|
| 86 | |
---|
| 87 | ! Duplicate on all longitudes: |
---|
| 88 | v1(2:, :, :) = spread(v1(1, :, :), dim=1, ncopies=iim-1) |
---|
| 89 | end if |
---|
| 90 | !$omp end master |
---|
| 91 | |
---|
| 92 | call scatter2d(v1, v2) |
---|
| 93 | |
---|
| 94 | ! Regrid in pressure at each horizontal position: |
---|
| 95 | do i = 1, klon |
---|
| 96 | v3(i, llm:1:-1) = regr1_step_av(v2(i, :), press_in_edg, & |
---|
| 97 | paprs(i, llm+1:1:-1)) |
---|
| 98 | ! (invert order of indices because "paprs" is in descending order) |
---|
| 99 | end do |
---|
| 100 | |
---|
| 101 | end subroutine regr_pr_av |
---|
| 102 | |
---|
| 103 | !*************************************************************** |
---|
| 104 | |
---|
| 105 | subroutine regr_pr_int(ncid, name, julien, plev, pplay, top_value, v3) |
---|
| 106 | |
---|
| 107 | ! "regr_pr_int" stands for "regrid pressure interpolation". |
---|
| 108 | ! The target vertical LMDZ grid is the grid of mid-layers. |
---|
| 109 | ! Regridding is by linear interpolation. |
---|
| 110 | |
---|
| 111 | use dimphy, only: klon |
---|
| 112 | use netcdf95, only: nf95_inq_varid, handle_err |
---|
| 113 | use netcdf, only: nf90_get_var |
---|
| 114 | use assert_m, only: assert |
---|
| 115 | use regr1_lint_m, only: regr1_lint |
---|
| 116 | use mod_phys_lmdz_mpi_data, only: is_mpi_root |
---|
| 117 | |
---|
| 118 | use mod_phys_lmdz_transfert_para, only: scatter2d |
---|
| 119 | ! (pack to the LMDZ horizontal "physics" grid and scatter) |
---|
| 120 | |
---|
| 121 | integer, intent(in):: ncid ! NetCDF ID of the file |
---|
| 122 | character(len=*), intent(in):: name ! of the NetCDF variable |
---|
| 123 | integer, intent(in):: julien ! jour julien, 1 <= julien <= 360 |
---|
| 124 | |
---|
| 125 | real, intent(in):: plev(:) |
---|
| 126 | ! (pressure level of input data, in Pa, in strictly ascending order) |
---|
| 127 | |
---|
| 128 | real, intent(in):: pplay(:, :) ! (klon, llm) |
---|
| 129 | ! (pression pour le mileu de chaque couche, en Pa) |
---|
| 130 | |
---|
| 131 | real, intent(in):: top_value |
---|
| 132 | ! (extra value of field at 0 pressure) |
---|
| 133 | |
---|
| 134 | real, intent(out):: v3(:, :) ! (klon, llm) |
---|
| 135 | ! (regridded field on the partial "physics" grid) |
---|
| 136 | ! ("v3(i, k)" is at longitude "xlon(i)", latitude |
---|
| 137 | ! "xlat(i)", middle of layer "k".) |
---|
| 138 | |
---|
| 139 | ! Variables local to the procedure: |
---|
| 140 | |
---|
| 141 | include "dimensions.h" |
---|
| 142 | integer varid, ncerr ! for NetCDF |
---|
| 143 | |
---|
| 144 | real v1(iim, jjm + 1, 0:size(plev)) |
---|
| 145 | ! (input field at day "julien", on the global "dynamics" horizontal grid) |
---|
| 146 | ! (First dimension is for longitude. |
---|
| 147 | ! The value is the same for all longitudes. |
---|
| 148 | ! "v1(:, j, k >=1)" is at latitude "rlatu(j)" and pressure "plev(k)".) |
---|
| 149 | |
---|
| 150 | real v2(klon, 0:size(plev)) |
---|
| 151 | ! (field scattered to the partial "physics" horizontal grid) |
---|
| 152 | ! "v2(i, k >= 1)" is at longitude "xlon(i)", latitude "xlat(i)" |
---|
| 153 | ! and pressure "plev(k)".) |
---|
| 154 | |
---|
| 155 | integer i |
---|
| 156 | |
---|
| 157 | !-------------------------------------------- |
---|
| 158 | |
---|
| 159 | call assert(shape(v3) == (/klon, llm/), "regr_pr_int v3") |
---|
| 160 | call assert(shape(pplay) == (/klon, llm/), "regr_pr_int pplay") |
---|
| 161 | |
---|
| 162 | !$omp master |
---|
| 163 | if (is_mpi_root) then |
---|
| 164 | call nf95_inq_varid(ncid, name, varid) |
---|
| 165 | |
---|
| 166 | ! Get data at the right day from the input file: |
---|
| 167 | ncerr = nf90_get_var(ncid, varid, v1(1, :, 1:), start=(/1, 1, julien/)) |
---|
| 168 | call handle_err("regr_pr_int nf90_get_var " // name, ncerr, ncid) |
---|
| 169 | ! Latitudes are in ascending order in the input file while |
---|
| 170 | ! "rlatu" is in descending order so we need to invert order: |
---|
| 171 | v1(1, :, 1:) = v1(1, jjm+1:1:-1, 1:) |
---|
| 172 | |
---|
| 173 | ! Complete "v1" with the value at 0 pressure: |
---|
| 174 | v1(1, :, 0) = top_value |
---|
| 175 | |
---|
| 176 | ! Duplicate on all longitudes: |
---|
| 177 | v1(2:, :, :) = spread(v1(1, :, :), dim=1, ncopies=iim-1) |
---|
| 178 | end if |
---|
| 179 | !$omp end master |
---|
| 180 | |
---|
| 181 | call scatter2d(v1, v2) |
---|
| 182 | |
---|
| 183 | ! Regrid in pressure at each horizontal position: |
---|
| 184 | do i = 1, klon |
---|
| 185 | v3(i, llm:1:-1) = regr1_lint(v2(i, :), (/0., plev/), pplay(i, llm:1:-1)) |
---|
| 186 | ! (invert order of indices because "pplay" is in descending order) |
---|
| 187 | end do |
---|
| 188 | |
---|
| 189 | end subroutine regr_pr_int |
---|
| 190 | |
---|
| 191 | end module regr_pr |
---|