source: LMDZ5/trunk/libf/phylmd/open_climoz_m.F90 @ 2791

Last change on this file since 2791 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: 3.2 KB
Line 
1! $Id$
2module open_climoz_m
3
4  implicit none
5
6contains
7
8  subroutine open_climoz(ncid, press_in_cen, press_in_edg, time_in)
9
10    ! This procedure should be called once per "gcm" run, by a single
11    ! thread of each MPI process.
12    ! The root MPI process opens "climoz_LMDZ.nc", reads the pressure
13    ! levels and the times and broadcasts them to the other processes.
14
15    ! We assume that, in "climoz_LMDZ.nc", the pressure levels are in hPa
16    ! and in strictly ascending order.
17
18    use netcdf95, only: nf95_open, nf95_close, nf95_gw_var, nf95_inq_varid
19    use netcdf, only: nf90_nowrite
20
21    use mod_phys_lmdz_mpi_data, only: is_mpi_root
22    use mod_phys_lmdz_mpi_transfert, only: bcast_mpi ! broadcast
23    USE phys_cal_mod, ONLY: calend, year_len, year_cur
24
25    ! OpenMP shares arguments:
26    integer, intent(out):: ncid ! of "climoz_LMDZ.nc", OpenMP shared
27
28    ! pressure levels for ozone climatology, in Pa, strictly ascending order
29    real, pointer:: press_in_cen(:) !--- at cells centers
30    real, pointer:: press_in_edg(:) !--- at the interfaces (pressure intervals)
31
32    real, pointer:: time_in(:)
33    ! times of ozone climatology records, in days since Jan. 1st 0h
34
35    ! Variables local to the procedure:
36
37    integer varid ! for NetCDF
38    integer n_plev ! number of pressure levels in the input data
39    integer n_time ! number of time records    in the input field
40    integer k
41    CHARACTER(LEN=80)  :: modname
42    CHARACTER(LEN=320) :: msg
43
44    !---------------------------------------
45
46    modname="open_climoz"
47    print *, "Call sequence information: "//TRIM(modname)
48
49    if (is_mpi_root) then
50       call nf95_open("climoz_LMDZ.nc", nf90_nowrite, ncid)
51
52       call nf95_inq_varid(ncid, "plev", varid)
53       call nf95_gw_var(ncid, varid, press_in_cen)
54       ! Convert from hPa to Pa because "paprs" and "pplay" are in Pa:
55       press_in_cen = press_in_cen * 100.
56       n_plev = size(press_in_cen)
57       CALL NF95_INQ_VARID(ncID, "time", varID)
58       CALL NF95_GW_VAR(ncid, varid, time_in)
59       n_time = SIZE(time_in)
60       IF(ALL([year_len,14]/=n_time)) THEN             !--- Check records number
61         WRITE(msg,'(a,3(i4,a))')TRIM(modname)//': expected input file could be&
62         & monthly or daily (14 or ',year_len,' records for year ',year_cur,' a&
63         &nd calendar "'//TRIM(calend)//'"), but ',n_time,' records were found'
64         CALL abort_physic(modname, msg, 1)
65       END IF
66
67    end if
68
69    call bcast_mpi(n_plev)
70    if (.not. is_mpi_root) allocate(press_in_cen(n_plev))
71    call bcast_mpi(press_in_cen)
72    call bcast_mpi(n_time)
73    if (.not. is_mpi_root) allocate(time_in(n_time))
74    call bcast_mpi(time_in)
75
76    ! Compute edges of pressure intervals:
77    allocate(press_in_edg(n_plev + 1))
78    if (is_mpi_root) then
79       press_in_edg(1) = 0.
80       ! We choose edges halfway in logarithm:
81       forall (k = 2:n_plev) press_in_edg(k) = sqrt(press_in_cen(k - 1) * press_in_cen(k))
82       press_in_edg(n_plev + 1) = huge(0.)
83       ! (infinity, but any value guaranteed to be greater than the
84       ! surface pressure would do)
85    end if
86    call bcast_mpi(press_in_edg)
87!    deallocate(press_in_cen) ! pointer
88
89  end subroutine open_climoz
90
91end module open_climoz_m
Note: See TracBrowser for help on using the repository browser.