source: LMDZ6/trunk/libf/obsolete/regr1_lint_m.F90 @ 5286

Last change on this file since 5286 was 2788, checked in by dcugnet, 8 years ago

Changes in ce0l about the way ozone forcing files are generated:

1) 3D raw input files "climoz.nc" are now handled.
2) Default behaviour is now to let the gcm interpolate in time online.

This helps to avoid huge forcing files (in particular for 3D fields).
In this case, the output files "climoz_LMDZ.nc" all have 14 records:

  • records 2-13 are obtained with records 1-12 of "climoz.nc".
  • records 1 and 14 are obtained respectively with:
    • record 12 of "climoz_m.nc" if available, of "climoz.nc" otherwise.
    • record 1 of "climoz_p.nc" if available, of "climoz.nc" otherwise.

3) If ok_daily_climoz key is TRUE, the time interpolation (one record

a day) is forced, using the 14 records described below.
This now depends on the calendar (it was on a 360 days basis only).

Changes in the gcm about the way zone forcing files are read/interpolated:

1) 3D horizontally interpolated "climoz_LMDZ.nc" files are now handled.
2) Daily files (already interpolated in time) are still handled, but their

number of records must match the expected number of days, that depends
on the calendar (records step is no longer 1/360 year).

3) 14 records monthly files are now handled (and prefered). This reduces

the I/O to a minimum and the aditional computational cost is low (simple
online linear time interpolation).

4) If adjust_tropopause key is TRUE, the input fields are stretched using

following method:

  • LMDZ dynamical tropopause is detected: Ptrop_lmdz = MAX ( P(Potential Vorticity==2PVU), P(theta==380K) )
  • file chemical tropopause is detected: Ptrop_file = P( tro3 == o3t ), where:

o3t = 91. + 28. * SIN(PI*(month-2)/6) (ppbV)

This formula comes from Thouret & al., ACP 6, 1033-1051, 2006.
The second term of the expression is multiplied by TANH(lat_deg/20.)
to account for latitude dependency.

  • File profile is streched in a +/- 5kms zone around the mean tropopause to ensure resulting tropopause matches the one of LMDZ. See procedure regr_pr_time_av for more details.
  • Property copyright set to
    Name of program: LMDZ
    Creation date: 1984
    Version: LMDZ5
    License: CeCILL version 2
    Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
    See the license file in the root directory
File size: 2.6 KB
Line 
1! $Id$
2module 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
19contains
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
98end module regr1_lint_m
Note: See TracBrowser for help on using the repository browser.