template.f90

How to create a new message from a template.

00001 ! Copyright 2005-2007 ECMWF
00002 ! 
00003 ! Licensed under the GNU Lesser General Public License which
00004 ! incorporates the terms and conditions of version 3 of the GNU
00005 ! General Public License.
00006 ! See LICENSE and gpl-3.0.txt for details.
00007 !
00008 !
00009 !  Description: how to create a new GRIB message from a template.
00010 !               
00011 !
00012 !
00013 !  Author: Anne Fouilloux
00014 !
00015 !
00016 !
00017 program template
00018   use grib_api
00019   implicit none  
00020   integer, parameter                            :: NbPressureLevel = 10
00021   integer                                       :: err
00022   integer                                       :: nx, ny, i
00023   integer                                       :: outfile
00024   integer                                       :: igrib
00025   integer                                       :: date, level
00026   real                                          :: missingValue
00027   integer                                       :: indicatorOfParameter
00028   double precision, dimension(:,:), allocatable :: field2D 
00029   integer, dimension(NbPressureLevel)           :: levels=(/1,2,10,50,70,100,200,500,850,1000/)
00030 
00031   ! Initialisation
00032 
00033   date = 20080104
00034   indicatorOfParameter = 130 ! temperature in table 128 used in the template
00035 
00036   !     a new grib message is loaded from an existing template
00037   !     templates are searched in a default template path (use grib_info
00038   !     to see where it is. The default template path can be changed by
00039   !     setting the environment variable GRIB_TEMPLATES_PATH
00040   call grib_new_from_template(igrib, "GRIB1")
00041 
00042   ! Open a new GRIB file to write new GRIB messages
00043   call grib_open_file(outfile, 'out.grib1','w')
00044 
00045   call grib_get(igrib,"numberOfPointsAlongAParallel", nx)
00046   
00047   call grib_get(igrib,"numberOfPointsAlongAMeridian",ny)
00048 
00049   call grib_get(igrib,'lev',level)
00050 
00051   allocate(field2D(nx,ny),stat=err)
00052 
00053   if (err .ne. 0) then
00054      print*, 'Failed to allocate ', nx*ny, ' values'
00055      STOP
00056   end if  
00057   call generate_field(field2D, missingValue)
00058  
00059 ! Same date and same parameters for all the levels
00060   call grib_set(igrib,"dataDate", date)
00061   call grib_set(igrib,"indicatorOfParameter", indicatorOfParameter)
00062 
00063   ! Write new GRIBs from the same message for different pressure levels 
00064   do i=1, size(levels)
00065      call grib_set(igrib,'lev',levels(i))
00066 
00067      ! compute the field for each pressure level     
00068      call generate_field(field2D, missingValue)
00069 
00070      ! set the value of missingValue
00071      call grib_set(igrib, 'missingValue',missingValue)
00072 
00073      call grib_set(igrib,"bitmapPresent",1)
00074 
00075      ! use pack to create 1D values
00076      call grib_set(igrib,'values',pack(field2D, mask=.true.))
00077  
00078   !     write the new message to a file
00079      call grib_write(igrib,outfile)
00080   end do
00081 
00082   call grib_release(igrib)
00083 
00084   call grib_close_file(outfile)
00085   deallocate(field2D)
00086 
00087 contains
00088 !======================================
00089 subroutine generate_field(gfield2D, missingValue)
00090  double precision, dimension(:,:) :: gfield2D
00091  real, intent(out)                :: missingValue
00092 
00093  integer, parameter               :: borders = 5
00094  integer                          :: i
00095  integer, dimension(2)            :: shape_gfield2D
00096 
00097  shape_gfield2D = shape(gfield2D)
00098 
00099  call random_number(gfield2D)
00100 ! define the value of missingValue
00101  missingValue = 99999
00102 ! just keep the boundaries for gfield2D
00103   gfield2D(borders+1:shape_gfield2D(1)-borders,borders+1:shape_gfield2D(2)-borders) = missingValue
00104 
00105 end subroutine generate_field
00106 !======================================
00107 
00108 end program template

Generated on Mon Sep 29 16:35:07 2008 for grib_api by  doxygen 1.5.4