trunk/libf/regr3_lint.f90

module regr3_lint_m

  ! Author: Lionel GUEZ

  implicit none

  interface regr3_lint
     ! Each procedure regrids by linear interpolation.
     ! The regridding operation is done on the third dimension of the
     ! input array.
     ! The difference betwwen the procedures is the rank of the first argument.
     module procedure regr33_lint, regr34_lint
  end interface

  private
  public regr3_lint

contains

  function regr33_lint(vs, xs, xt) result(vt)

    ! "vs" has rank 3.

    use nr_util, only: assert_eq
    use numer_rec, only: hunt

    real, intent(in):: vs(:, :, :)
    ! (values of the function at source points "xs")

    real, intent(in):: xs(:)
    ! (abscissas of points in source grid, in strictly monotonic order)

    real, intent(in):: xt(:)
    ! (abscissas of points in target grid)

    real vt(size(vs, 1), size(vs, 2), size(xt))
    ! (values of the function on the target grid)

    ! Variables local to the procedure:
    integer is, it, ns
    integer is_b ! "is" bound between 1 and "ns - 1"

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

    ns = assert_eq(size(vs, 3), size(xs), "regr33_lint ns")

    is = -1 ! go immediately to bisection on first call to "hunt"

    do it = 1, size(xt)
       call hunt(xs, xt(it), is)
       is_b = min(max(is, 1), ns - 1)
       vt(:, :, it) = ((xs(is_b+1) - xt(it)) * vs(:, :, is_b) &
            + (xt(it) - xs(is_b)) * vs(:, :, is_b+1)) / (xs(is_b+1) - xs(is_b))
    end do

  end function regr33_lint

  !*********************************************************

  function regr34_lint(vs, xs, xt) result(vt)

    ! "vs" has rank 4.

    use nr_util, only: assert_eq
    use numer_rec, only: hunt

    real, intent(in):: vs(:, :, :, :)
    ! (values of the function at source points "xs")

    real, intent(in):: xs(:)
    ! (abscissas of points in source grid, in strictly monotonic order)

    real, intent(in):: xt(:)
    ! (abscissas of points in target grid)

    real vt(size(vs, 1), size(vs, 2), size(xt), size(vs, 4))
    ! (values of the function on the target grid)

    ! Variables local to the procedure:
    integer is, it, ns
    integer is_b ! "is" bound between 1 and "ns - 1"

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

    ns = assert_eq(size(vs, 3), size(xs), "regr34_lint ns")

    is = -1 ! go immediately to bisection on first call to "hunt"

    do it = 1, size(xt)
       call hunt(xs, xt(it), is)
       is_b = min(max(is, 1), ns - 1)
       vt(:, :, it, :) = ((xs(is_b+1) - xt(it)) * vs(:, :, is_b, :) &
            + (xt(it) - xs(is_b)) * vs(:, :, is_b+1, :)) &
            / (xs(is_b+1) - xs(is_b))
    end do

  end function regr34_lint

end module regr3_lint_m
1	guez	7	module regr3_lint_m
2
3	guez	36	! Author: Lionel GUEZ
4
5	guez	7	implicit none
6
7			interface regr3_lint
8			! Each procedure regrids by linear interpolation.
9			! The regridding operation is done on the third dimension of the
10			! input array.
11			! The difference betwwen the procedures is the rank of the first argument.
12	guez	36	module procedure regr33_lint, regr34_lint
13	guez	7	end interface
14
15			private
16			public regr3_lint
17
18			contains
19
20			function regr33_lint(vs, xs, xt) result(vt)
21
22			! "vs" has rank 3.
23
24	guez	36	use nr_util, only: assert_eq
25			use numer_rec, only: hunt
26	guez	7
27			real, intent(in):: vs(:, :, :)
28			! (values of the function at source points "xs")
29
30			real, intent(in):: xs(:)
31			! (abscissas of points in source grid, in strictly monotonic order)
32
33			real, intent(in):: xt(:)
34			! (abscissas of points in target grid)
35
36			real vt(size(vs, 1), size(vs, 2), size(xt))
37			! (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, 3), size(xs), "regr33_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			vt(:, :, it) = ((xs(is_b+1) - xt(it)) * vs(:, :, is_b) &
53			+ (xt(it) - xs(is_b)) * vs(:, :, is_b+1)) / (xs(is_b+1) - xs(is_b))
54			end do
55
56			end function regr33_lint
57
58	guez	36	!*********************************************************
59
60			function regr34_lint(vs, xs, xt) result(vt)
61
62			! "vs" has rank 4.
63
64			use nr_util, only: assert_eq
65			use numer_rec, only: hunt
66
67			real, intent(in):: vs(:, :, :, :)
68			! (values of the function at source points "xs")
69
70			real, intent(in):: xs(:)
71			! (abscissas of points in source grid, in strictly monotonic order)
72
73			real, intent(in):: xt(:)
74			! (abscissas of points in target grid)
75
76			real vt(size(vs, 1), size(vs, 2), size(xt), size(vs, 4))
77			! (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, 3), size(xs), "regr34_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, :)) &
94			/ (xs(is_b+1) - xs(is_b))
95			end do
96
97			end function regr34_lint
98
99	guez	7	end module regr3_lint_m