1! Copyright 2005-2018 ECMWF.
2!
3! This software is licensed under the terms of the Apache Licence Version 2.0
4! which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
5!
6! In applying this licence, ECMWF does not waive the privileges and immunities granted to it by
7! virtue of its status as an intergovernmental organisation nor does it submit to any jurisdiction.
8!
9!
10!
11!  Description: how to set pv values.
12!
13!
14program set_pv
15  use grib_api
16  implicit none
17  integer                         :: numberOfLevels
18  integer                         :: numberOfCoefficients
19  integer                         :: outfile, igrib
20  integer                         :: i, ios
21  real, dimension(:),allocatable  :: pv
22
23  numberOfLevels=60
24  numberOfCoefficients=2*(numberOfLevels+1)
25
26  allocate(pv(numberOfCoefficients))
27
28  ! read the model level coefficients from file
29  open( unit=1, file="../../data/60_model_levels", &
30                form="formatted",action="read")
31
32  do i=1,numberOfCoefficients,2
33    read(unit=1,fmt=*, iostat=ios) pv(i), pv(i+1)
34    if (ios /= 0) then
35      print *, "I/O error: ",ios
36      exit
37    end if
38  end do
39
40  ! print coefficients
41  !do i=1,numberOfCoefficients,2
42  !  print *,"  a=",pv(i)," b=",pv(i+1)
43  !end do
44
45  close(unit=1)
46
47  call grib_open_file(outfile, 'out.pv.grib1','w')
48
49  !     a new grib message is loaded from file
50  !     igrib is the grib id to be used in subsequent calls
51  call grib_new_from_samples(igrib, "reduced_gg_sfc_grib1")
52
53  !     set levtype to ml (model level)
54  call grib_set(igrib,'typeOfLevel','hybrid')
55
56  !     set level
57  call grib_set(igrib,'level',2)
58
59  !     set PVPresent as an integer
60  call grib_set(igrib,'PVPresent',1)
61
62  call grib_set(igrib,'pv',pv)
63
64  !     write modified message to a file
65  call grib_write(igrib,outfile)
66
67  !  FREE MEMORY
68  call grib_release(igrib)
69  deallocate(pv)
70
71  call grib_close_file(outfile)
72
73end program set_pv
74