samples.f90

How to create a new message from a samples.

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 sample.
00010 !               
00011 !
00012 program sample
00013   use grib_api
00014   implicit none  
00015   integer  :: err
00016   integer  :: outfile, infile, datafile
00017   integer  :: igribsample,igribclone,igribdata, size
00018   integer  :: date, startStep, endStep, table2Version, indicatorOfParameter
00019   integer  :: decimalPrecision
00020   real     :: missingValue
00021   character(len=10) stepType
00022   double precision, dimension(:), allocatable   :: v1,v2,v 
00023 
00024   date = 20080104
00025   startStep = 0
00026   endStep = 12
00027   stepType = 'accum'
00028   table2Version = 2
00029   indicatorOfParameter = 61
00030   decimalPrecision = 2
00031 
00032   !     a new grib message is loaded from an existing sample
00033   !     samples are searched in a default sample path (use grib_info
00034   !     to see where it is. The default sample path can be changed by
00035   !     setting the environment variable GRIB_SAMPLES_PATH
00036   call grib_new_from_samples(igribsample, "regular_latlon_surface.grib1")
00037 
00038   call grib_open_file(outfile, 'out.grib1','w')
00039   call grib_open_file(datafile,'../../data/tp_ecmwf.grib','r')
00040 
00041   call grib_new_from_file(datafile,igribdata,err)
00042 
00043   call grib_get_size(igribdata,'values',size)
00044   allocate(v(size))
00045   allocate(v1(size))
00046   allocate(v2(size))
00047 
00048   call grib_get(igribdata,'values',v)
00049 
00050   v=v*1000.0 ! different units for the output grib
00051   v1=v
00052 
00053   do while (err/=GRIB_END_OF_FILE) 
00054  
00055     call grib_clone(igribsample,igribclone) ! clone sample before modifying it
00056 
00057     call grib_set(igribclone,'date',date)
00058     call grib_set(igribclone,'table2Version',table2Version)
00059     call grib_set(igribclone,'indicatorOfParameter',indicatorOfParameter)
00060 
00061     call grib_set(igribclone,'stepType',stepType)
00062     call grib_set(igribclone,'startStep',startStep)
00063     call grib_set(igribclone,'endStep',endStep)
00064 
00065     call grib_set(igribclone,'decimalPrecision',decimalPrecision)
00066 
00067     call grib_set(igribclone,'values',v)
00068 
00069     call grib_write(igribclone,outfile)
00070     
00071     call grib_new_from_file(datafile,igribdata,err)
00072 
00073     if (err==0) then
00074       call grib_get(igribdata,'values',v2)
00075 
00076       v2=v2*1000.0 ! different units for the output grib
00077 
00078       v=v2-v1 ! accumulation from startStep to endStep
00079 
00080       v1=v2 ! save previous step field
00081 
00082       startStep=startStep+12
00083       endStep=endStep+12
00084 
00085     endif
00086 
00087   enddo
00088 
00089   deallocate(v)
00090   deallocate(v1)
00091   deallocate(v2)
00092 
00093   call grib_close_file(outfile)
00094 
00095 end program sample

Generated on Tue Sep 22 15:18:22 2009 for grib_api by  doxygen 1.5.3