! $Id$

module regr_lat_time_coefoz_m

  ! Author: Lionel GUEZ

  implicit none

  private

  public regr_lat_time_coefoz

contains

    ! "regr_lat_time_coefoz" stands for "regrid latitude time

    ! coefficients ozone".

    ! This procedure reads from a NetCDF file coefficients for ozone

    ! chemistry, regrids them in latitude and time, and writes the

    ! regridded fields to a new NetCDF file.

    ! The input fields depend on time, pressure level and latitude.

    ! We assume that the input fields are step functions of latitude.

    ! Regridding in latitude is made by averaging, with a cosine of

    ! latitude factor.

    ! The target LMDZ latitude grid is the "scalar" grid: "rlatu".

    ! The values of "rlatu" are taken to be the centers of intervals.

    ! Regridding in time is by linear interpolation.

    ! Monthly values are processed to get daily values, on the basis

    ! of a 360-day calendar.

    ! We assume that in the input file:

    ! -- the latitude is in degrees and strictly monotonic (as all

    ! NetCDF coordinate variables should be);

    ! -- time increases from January to December (even though we do

    ! not use values of the input time coordinate);

    ! -- pressure is in hPa and in strictly ascending order (even

    ! though we do not use pressure values here, we write the unit of

    ! pressure in the NetCDF header, and we will use the assumptions later,

    ! when we regrid in pressure).

    use mod_grid_phy_lmdz, ONLY : nbp_lat

    use regr_conserv_m, only: regr_conserv

    use regr_lint_m, only: regr_lint

    use netcdf95, only: nf95_open, nf95_close, nf95_inq_varid, nf95_get_var, &

         nf95_put_var, nf95_gw_var

    use netcdf, only: nf90_nowrite

    use nrtype, only: pi

    use regular_lonlat_mod, only: boundslat_reg, south

    ! Variables local to the procedure:

    integer ncid_in, ncid_out ! NetCDF IDs for input and output files

    integer n_plev ! number of pressure levels in the input data

    integer n_lat! number of latitudes in the input data

    ! (of input data, converted to rad, sorted in strictly ascending order)

    real, allocatable:: lat_in_edg(:)

    ! (edges of latitude intervals for input data, in rad, in strictly

    ! ascending order)

    logical desc_lat ! latitude in descending order in the input file

    ! (ozone parameter from the input file)

    ! ("o3_par_in(j, l, month)" is at latitude "latitude(j)" and pressure

    ! level "plev(l)". "month" is between 1 and 12.)

    ! (mean of a variable "v" over a latitude interval)

    ! (First dimension is latitude interval.

    ! The latitude interval for "v_regr_lat(j,:, :)" contains "rlatu(j)".

    ! If "j" is between 2 and "jjm" then the interval is:

    ! [rlatv(j), rlatv(j-1)]

    ! If "j" is 1 or "jjm + 1" then the interval is:

    ! [rlatv(1), pi / 2]

    ! or:

    ! [- pi / 2, rlatv(jjm)]

    ! respectively.

    ! "v_regr_lat(:, l, :)" is for pressure level "plev(l)".

    ! Last dimension is month number.)

    ! (regridded ozone parameter)

    ! ("o3_par_out(j, l, day)" is at latitude "rlatu(j)", pressure

    ! level "plev(l)" and date "January 1st 0h" + "tmidday(day)", in a

    ! 360-day calendar.)

    integer j

    integer i_v ! index of ozone parameter

    integer, parameter:: n_o3_param = 8 ! number of ozone parameters

    character(len=11) name_in(n_o3_param)

    ! (name of NetCDF primary variable in the input file)

    character(len=9) name_out(n_o3_param)

    ! (name of NetCDF primary variable in the output file)

    integer varid_in(n_o3_param), varid_out(n_o3_param), varid_plev, varid_time

    integer ncerr, varid

    ! (for NetCDF)

    real, parameter:: tmidmonth(0:13) = (/(-15. + 30. * j, j = 0, 13)/)

    ! (time to middle of month, in days since January 1st 0h, in a

    ! 360-day calendar)

    ! (We add values -15 and 375 so that, for example, day 3 of the year is

    ! interpolated between the December and the January value.)

    real, parameter:: tmidday(360) = (/(j + 0.5, j = 0, 359)/)

    ! (time to middle of day, in days since January 1st 0h, in a

    ! 360-day calendar)

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

    ! Names of ozone parameters:

    i_v = 0

    i_v = i_v + 1

    i_v = i_v + 1

    i_v = i_v + 1

    i_v = i_v + 1

    i_v = i_v + 1

    i_v = i_v + 1

    i_v = i_v + 1

    i_v = i_v + 1

    ! Get coordinates from the input file:

    call nf95_gw_var(ncid_in, varid, latitude)

    ! Convert from degrees to rad, because "boundslat_reg" is in rad:

    ! We need to supply the latitudes to "regr_conserv" in

    ! ascending order, so invert order if necessary:

    ! Compute edges of latitude intervals:

    call nf95_gw_var(ncid_in, varid, plev)

    ! (We only need the pressure coordinate to copy it to the output file.)

    ! Get the IDs of ozone parameters in the input file:

    end do

    ! Create the output file and get the variable IDs:

    call prepare_out(ncid_in, varid_in, n_plev, name_out, ncid_out, &

    ! Write remaining coordinate variables:

       ! Process ozone parameter "name_in(i_v)"

       ! Regrid in latitude:

       ! We average with respect to sine of latitude, which is

       ! equivalent to weighting by cosine of latitude:

       call regr_conserv(1, o3_par_in, xs = sin(lat_in_edg), &

            xt = (/-1., sin((/boundslat_reg(nbp_lat-1:1:-1,south)/)), 1./), &

       ! (invert order of indices in "v_regr_lat" because "rlatu" is

       ! in descending order)

       ! Duplicate January and December values, in preparation of time

       ! interpolation:

       ! Regrid in time by linear interpolation:

       ! Write to file:

       call nf95_put_var(ncid_out, varid_out(i_v), &

       ! (The order of "rlatu" is inverted in the output file)

    end do

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

    ! This subroutine creates the NetCDF output file, defines

    ! dimensions and variables, and writes one of the coordinate variables.

    use mod_grid_phy_lmdz, ONLY : nbp_lat

    use assert_eq_m, only: assert_eq

    use netcdf95, only: nf95_create, nf95_def_dim, nf95_def_var, &

         nf95_put_att, nf95_enddef, nf95_copy_att, nf95_put_var

    use netcdf, only: nf90_clobber, nf90_float, nf90_copy_att, nf90_global

    use nrtype, only: pi

    use regular_lonlat_mod, only : lat_reg

    integer, intent(in):: ncid_in, varid_in(:), n_plev

    character(len=*), intent(in):: name_out(:) ! of NetCDF variables

    integer, intent(out):: ncid_out, varid_out(:), varid_plev, varid_time

    ! Variables local to the procedure:

    integer ncerr

    integer dimid_rlatu, dimid_plev, dimid_time

    integer varid_rlatu

    integer i, n_o3_param

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

    n_o3_param = assert_eq(size(varid_in), size(varid_out), &

    ! Dimensions:

    ! Define coordinate variables:

    ! Define primary variables:

       call nf95_def_var(ncid_out, name_out(i), nf90_float, &

       ! The following commands may fail. That is OK. It should just

       ! mean that the attribute is not defined in the input file.

       ncerr = nf90_copy_att(ncid_in, varid_in(i), "long_name",&

       ncerr = nf90_copy_att(ncid_in, varid_in(i), "units", ncid_out,&

       ncerr = nf90_copy_att(ncid_in, varid_in(i), "standard_name", ncid_out,&

    end do

    ! Global attributes:

    call nf95_copy_att(ncid_in, nf90_global, "Conventions", ncid_out, &

    ! Write one of the coordinate variables:

    ! (convert from rad to degrees and sort in ascending order)

  contains

      use netcdf, only: nf90_noerr, nf90_strerror

      character(len=*), intent(in):: att_name

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

         print *, "prepare_out " // trim(name_out(i)) &

              // " nf90_copy_att " // att_name // " -- " &

              // trim(nf90_strerror(ncerr))

      end if

  end subroutine prepare_out

end module regr_lat_time_coefoz_m