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: count messages before processing 00010 ! 00011 ! Author: Enrico Fucile 00012 ! 00013 ! 00014 program get 00015 use grib_api 00016 implicit none 00017 00018 integer :: ifile 00019 integer :: iret 00020 integer :: n 00021 integer :: i 00022 integer,dimension(:),allocatable :: igrib 00023 real :: latitudeOfFirstPointInDegrees 00024 real :: longitudeOfFirstPointInDegrees 00025 real :: latitudeOfLastPointInDegrees 00026 real :: longitudeOfLastPointInDegrees 00027 integer :: numberOfPointsAlongAParallel 00028 integer :: numberOfPointsAlongAMeridian 00029 real, dimension(:), allocatable :: values 00030 integer :: numberOfValues 00031 real :: average,min_val, max_val 00032 00033 call grib_open_file(ifile, & 00034 '../../data/tigge_pf_ecmwf.grib2','r') 00035 00036 ! count the messages in the file 00037 call grib_count_in_file(ifile,n) 00038 allocate(igrib(n)) 00039 igrib=-1 00040 00041 ! Load the messages from the file. 00042 DO i=1,n 00043 call grib_new_from_file(ifile,igrib(i), iret) 00044 END DO 00045 00046 ! we can close the file 00047 call grib_close_file(ifile) 00048 00049 ! Loop on all the messages in memory 00050 DO i=1,n 00051 write(*,*) 'processing message number ',i 00052 ! get as a integer 00053 call grib_get(igrib(i),'numberOfPointsAlongAParallel', & 00054 numberOfPointsAlongAParallel) 00055 write(*,*) 'numberOfPointsAlongAParallel=', & 00056 numberOfPointsAlongAParallel 00057 00058 ! get as a integer 00059 call grib_get(igrib(i),'numberOfPointsAlongAMeridian', & 00060 numberOfPointsAlongAMeridian) 00061 write(*,*) 'numberOfPointsAlongAMeridian=', & 00062 numberOfPointsAlongAMeridian 00063 00064 ! get as a real 00065 call grib_get(igrib(i), 'latitudeOfFirstGridPointInDegrees', & 00066 latitudeOfFirstPointInDegrees) 00067 write(*,*) 'latitudeOfFirstGridPointInDegrees=', & 00068 latitudeOfFirstPointInDegrees 00069 00070 ! get as a real 00071 call grib_get(igrib(i), 'longitudeOfFirstGridPointInDegrees', & 00072 longitudeOfFirstPointInDegrees) 00073 write(*,*) 'longitudeOfFirstGridPointInDegrees=', & 00074 longitudeOfFirstPointInDegrees 00075 00076 ! get as a real 00077 call grib_get(igrib(i), 'latitudeOfLastGridPointInDegrees', & 00078 latitudeOfLastPointInDegrees) 00079 write(*,*) 'latitudeOfLastGridPointInDegrees=', & 00080 latitudeOfLastPointInDegrees 00081 00082 ! get as a real 00083 call grib_get(igrib(i), 'longitudeOfLastGridPointInDegrees', & 00084 longitudeOfLastPointInDegrees) 00085 write(*,*) 'longitudeOfLastGridPointInDegrees=', & 00086 longitudeOfLastPointInDegrees 00087 00088 00089 ! get the size of the values array 00090 call grib_get_size(igrib(i),'values',numberOfValues) 00091 write(*,*) 'numberOfValues=',numberOfValues 00092 00093 allocate(values(numberOfValues), stat=iret) 00094 ! get data values 00095 call grib_get(igrib(i),'values',values) 00096 call grib_get(igrib(i),'min',min_val) ! can also be obtained through minval(values) 00097 call grib_get(igrib(i),'max',max_val) ! can also be obtained through maxval(values) 00098 call grib_get(igrib(i),'average',average) ! can also be obtained through maxval(values) 00099 00100 write(*,*)'There are ',numberOfValues, & 00101 ' average is ',average, & 00102 ' min is ', min_val, & 00103 ' max is ', max_val 00104 write(*,*) '---------------------' 00105 END DO 00106 00107 DO i=1,n 00108 call grib_release(igrib(i)) 00109 END DO 00110 00111 deallocate(values) 00112 deallocate(igrib) 00113 00114 end program get