4 |
|
|
5 |
contains |
contains |
6 |
|
|
7 |
subroutine inifilr_hemisph(rlat, colat0, rlamda, unit, eignfn, jfilt, & |
subroutine inifilr_hemisph(rlat, rlamda, unit, eignfn, jfilt, matrice, & |
8 |
matrice, matrinv) |
matrinv) |
9 |
|
|
10 |
! See notes. |
! See notes. |
11 |
|
|
12 |
USE dimens_m, ONLY : iim |
USE dimensions, ONLY : iim |
13 |
use nr_util, only: pi, ifirstloc |
use nr_util, only: pi, ifirstloc |
14 |
|
|
15 |
real, intent(in):: rlat(:) ! (n) |
real, intent(in):: rlat(:) ! (n_lat) |
16 |
! latitudes, in rad, in [0, pi / 2[, in strictly ascending order |
! latitudes, in rad, in [0, pi / 2[, in strictly ascending order |
17 |
|
|
|
REAL, intent(in):: colat0 ! > 0 |
|
18 |
REAL, intent(in):: rlamda(2:) ! (2:iim) > 0, in descending order |
REAL, intent(in):: rlamda(2:) ! (2:iim) > 0, in descending order |
19 |
integer, intent(in):: unit |
integer, intent(in):: unit |
20 |
|
|
22 |
! eigenvectors of the discrete second derivative with respect to longitude |
! eigenvectors of the discrete second derivative with respect to longitude |
23 |
|
|
24 |
integer, intent(out):: jfilt |
integer, intent(out):: jfilt |
|
real, pointer:: matrice(:, :, :) ! (iim, iim, n - jfilt + 1) matrice filtre |
|
25 |
|
|
26 |
real, pointer, optional:: matrinv(:, :, :) ! (iim, iim, n - jfilt + 1) |
real, pointer:: matrice(:, :, :) ! (iim, iim, n_lat - jfilt + 1) |
27 |
|
! matrice filtre |
28 |
|
|
29 |
|
real, pointer, optional:: matrinv(:, :, :) ! (iim, iim, n_lat - jfilt + 1) |
30 |
! matrice pour le filtre inverse |
! matrice pour le filtre inverse |
31 |
|
|
32 |
! Local: |
! Local: |
33 |
|
|
34 |
integer n, i, j |
integer n_lat, i, j |
35 |
REAL eignft(iim, iim) |
REAL eignft(iim, iim) |
36 |
|
|
37 |
! Index of the mode from where modes are filtered: |
! Index of the mode from where modes are filtered: |
38 |
integer, allocatable:: modfrst(:) ! (jfilt:n) in {2, ..., iim} |
integer modfrst ! in {2, ..., iim} |
39 |
|
|
40 |
! Filtering coefficients (lamda_max * cos(rlat) / lamda): |
! Filtering coefficients (lamda_max * cos(rlat) / lamda): |
41 |
real coefil(2:iim) |
real coefil(2:iim) |
42 |
|
|
43 |
!----------------------------------------------------------- |
!----------------------------------------------------------- |
44 |
|
|
45 |
n = size(rlat) |
n_lat = size(rlat) |
46 |
jfilt = ifirstloc(cos(rlat) < min(colat0, 1. / rlamda(iim))) |
jfilt = ifirstloc(cos(rlat) < 1. / rlamda(iim)) |
47 |
allocate(modfrst(jfilt:n)) |
allocate(matrice(iim, iim, n_lat - jfilt + 1)) |
48 |
|
if (present(matrinv)) allocate(matrinv(iim, iim, n_lat - jfilt + 1)) |
49 |
DO j = jfilt, n |
|
50 |
modfrst(j) = ifirstloc(rlamda < 1. / cos(rlat(j)), my_lbound = 2) |
DO j = jfilt, n_lat |
51 |
write(unit, fmt = *) rlat(j) / pi * 180., modfrst(j) |
modfrst = ifirstloc(rlamda < 1. / cos(rlat(j)), my_lbound = 2) |
52 |
END DO |
write(unit, fmt = *) rlat(j) / pi * 180., modfrst |
53 |
|
coefil(modfrst:) = rlamda(modfrst:) * cos(rlat(j)) - 1. |
54 |
allocate(matrice(iim, iim, n - jfilt + 1)) |
eignft(:modfrst - 1, :) = 0. |
|
if (present(matrinv)) allocate(matrinv(iim, iim, n - jfilt + 1)) |
|
|
|
|
|
DO j = jfilt, n |
|
|
coefil(modfrst(j):) = rlamda(modfrst(j):) * cos(rlat(j)) - 1. |
|
|
eignft(:modfrst(j) - 1, :) = 0. |
|
55 |
|
|
56 |
forall (i = modfrst(j):iim) eignft(i, :) = eignfn(:, i) * coefil(i) |
forall (i = modfrst:iim) eignft(i, :) = eignfn(:, i) * coefil(i) |
57 |
matrice(:, :, j - jfilt + 1) = matmul(eignfn, eignft) |
matrice(:, :, j - jfilt + 1) = matmul(eignfn, eignft) |
58 |
|
|
59 |
if (present(matrinv)) then |
if (present(matrinv)) then |
60 |
forall (i = modfrst(j):iim) eignft(i, :) = eignfn(:, i) * coefil(i) & |
forall (i = modfrst:iim) eignft(i, :) = eignfn(:, i) * coefil(i) & |
61 |
/ (1. + coefil(i)) |
/ (1. + coefil(i)) |
62 |
matrinv(:, :, j - jfilt + 1) = matmul(eignfn, eignft) |
matrinv(:, :, j - jfilt + 1) = matmul(eignfn, eignft) |
63 |
end if |
end if |