1 | ! $Id$ |
---|
2 | module regr1_lint_m |
---|
3 | |
---|
4 | ! Author: Lionel GUEZ |
---|
5 | |
---|
6 | implicit none |
---|
7 | |
---|
8 | interface regr1_lint |
---|
9 | ! Each procedure regrids by linear interpolation. |
---|
10 | ! The regridding operation is done on the first dimension of the |
---|
11 | ! input array. |
---|
12 | ! The difference betwwen the procedures is the rank of the first argument. |
---|
13 | module procedure regr11_lint, regr12_lint |
---|
14 | end interface |
---|
15 | |
---|
16 | private |
---|
17 | public regr1_lint |
---|
18 | |
---|
19 | contains |
---|
20 | |
---|
21 | function regr11_lint(vs, xs, xt) result(vt) |
---|
22 | |
---|
23 | ! "vs" has rank 1. |
---|
24 | |
---|
25 | use assert_eq_m, only: assert_eq |
---|
26 | use interpolation, only: hunt !!, polint |
---|
27 | |
---|
28 | real, intent(in):: vs(:) |
---|
29 | ! (values of the function at source points "xs") |
---|
30 | |
---|
31 | real, intent(in):: xs(:) |
---|
32 | ! (abscissas of points in source grid, in strictly monotonic order) |
---|
33 | |
---|
34 | real, intent(in):: xt(:) |
---|
35 | ! (abscissas of points in target grid) |
---|
36 | |
---|
37 | real vt(size(xt)) ! values of the function on the target grid |
---|
38 | |
---|
39 | ! Variables local to the procedure: |
---|
40 | integer is, it, ns |
---|
41 | integer is_b ! "is" bound between 1 and "ns - 1" |
---|
42 | |
---|
43 | !-------------------------------------- |
---|
44 | |
---|
45 | ns = assert_eq(size(vs), size(xs), "regr11_lint ns") |
---|
46 | |
---|
47 | is = -1 ! go immediately to bisection on first call to "hunt" |
---|
48 | |
---|
49 | do it = 1, size(xt) |
---|
50 | call hunt(xs, xt(it), is) |
---|
51 | is_b = min(max(is, 1), ns - 1) |
---|
52 | !! call polint(xs(is_b:is_b+1), vs(is_b:is_b+1), xt(it), vt(it)) |
---|
53 | vt(it) = ((xs(is_b+1) - xt(it)) * vs(is_b) & |
---|
54 | + (xt(it) - xs(is_b)) * vs(is_b+1)) / (xs(is_b+1) - xs(is_b)) |
---|
55 | END DO |
---|
56 | |
---|
57 | end function regr11_lint |
---|
58 | |
---|
59 | !********************************************************* |
---|
60 | |
---|
61 | function regr12_lint(vs, xs, xt) result(vt) |
---|
62 | |
---|
63 | ! "vs" has rank 2. |
---|
64 | |
---|
65 | use assert_eq_m, only: assert_eq |
---|
66 | use interpolation, only: hunt |
---|
67 | |
---|
68 | real, intent(in):: vs(:, :) |
---|
69 | ! (values of the function at source points "xs") |
---|
70 | |
---|
71 | real, intent(in):: xs(:) |
---|
72 | ! (abscissas of points in source grid, in strictly monotonic order) |
---|
73 | |
---|
74 | real, intent(in):: xt(:) |
---|
75 | ! (abscissas of points in target grid) |
---|
76 | |
---|
77 | real vt(size(xt), size(vs, 2)) ! values of the function on the target grid |
---|
78 | |
---|
79 | ! Variables local to the procedure: |
---|
80 | integer is, it, ns |
---|
81 | integer is_b ! "is" bound between 1 and "ns - 1" |
---|
82 | |
---|
83 | !-------------------------------------- |
---|
84 | |
---|
85 | ns = assert_eq(size(vs, 1), size(xs), "regr12_lint ns") |
---|
86 | |
---|
87 | is = -1 ! go immediately to bisection on first call to "hunt" |
---|
88 | |
---|
89 | do it = 1, size(xt) |
---|
90 | call hunt(xs, xt(it), is) |
---|
91 | is_b = min(max(is, 1), ns - 1) |
---|
92 | vt(it, :) = ((xs(is_b+1) - xt(it)) * vs(is_b, :) & |
---|
93 | + (xt(it) - xs(is_b)) * vs(is_b+1, :)) / (xs(is_b+1) - xs(is_b)) |
---|
94 | END DO |
---|
95 | |
---|
96 | end function regr12_lint |
---|
97 | |
---|
98 | end module regr1_lint_m |
---|