From cfef150e1952998c5b087578fe57c2107ea61db5 Mon Sep 17 00:00:00 2001 From: Shahram Najm Date: Fri, 11 Dec 2015 14:14:41 +0000 Subject: [PATCH] ECC-196: Update the bufr_read_temp example --- data/bufr/bufr_data_files.txt | 1 + examples/F90/bufr_read_temp.f90 | 132 +++++++++++++----------------- examples/F90/bufr_read_temp.sh | 9 +- examples/python/bufr_read_temp.py | 101 +++++++---------------- examples/python/bufr_read_temp.sh | 10 +-- 5 files changed, 89 insertions(+), 164 deletions(-) diff --git a/data/bufr/bufr_data_files.txt b/data/bufr/bufr_data_files.txt index 7d8d0df21..e05afb09b 100644 --- a/data/bufr/bufr_data_files.txt +++ b/data/bufr/bufr_data_files.txt @@ -136,3 +136,4 @@ wavb_134.bufr new.bufr tropical_cyclone.bufr uegabe.bufr +PraticaTemp.bufr diff --git a/examples/F90/bufr_read_temp.f90 b/examples/F90/bufr_read_temp.f90 index 7134129c6..5fc8fd2e0 100644 --- a/examples/F90/bufr_read_temp.f90 +++ b/examples/F90/bufr_read_temp.f90 @@ -1,103 +1,81 @@ -! -!Copyright 2005-2015 ECMWF. +! Copyright 2005-2015 ECMWF. ! ! This software is licensed under the terms of the Apache Licence Version 2.0 -!which can be obtained at http://www.apache.org/licenses/LICENSE-2.0. +! which can be obtained at http://www.apache.org/licenses/LICENSE-2.0. ! ! In applying this licence, ECMWF does not waive the privileges and immunities granted to it by ! virtue of its status as an intergovernmental organisation nor does it submit to any jurisdiction. ! ! -! FOTRAN 90 Implementation: bufr_read_temp +! FORTRAN 90 Implementation: bufr_read_temp ! -! Description: how to read temperature significant levels from TEMP BUFR messages. +! Description: how to read levels from TEMP BUFR messages. ! ! Please note that TEMP reports can be encoded in various ways in BUFR. Therefore the code ! below might not work directly for other types of TEMP messages than the one used in the ! example. It is advised to use bufr_dump first to understand the structure of these messages. - +! program bufr_read_temp -use eccodes -implicit none -integer :: ifile -integer :: iret -integer :: ibufr -integer :: i, count=0 -integer(kind=4) :: num -real(kind=8), dimension(:), allocatable :: presVal,geoVal,tVal ,tdVal -character(len=128) :: keyName + use eccodes + implicit none + integer :: ifile + integer :: iret,ierr + integer :: ibufr + integer :: i, count=0 + integer(kind=4),dimension(:), allocatable :: timePeriod,extendedVerticalSoundingSignificance + integer(kind=4) :: blockNumber,stationNumber,numberOfLevels + character(len=30) :: str + real(kind=8),dimension(:), allocatable :: pressure,airTemperature,dewpointTemperature + real(kind=8),dimension(:), allocatable :: geopotentialHeight,latitudeDisplacement,longitudeDisplacement + real(kind=8),dimension(:), allocatable :: windDirection,windSpeed + character(len=128) :: keyName - call codes_open_file(ifile,'../../data/bufr/temp_101.bufr','r') + call codes_open_file(ifile,'../../data/bufr/PraticaTemp.bufr','r') -! the first bufr message is loaded from file -! ibufr is the bufr id to be used in subsequent calls + ! the first bufr message is loaded from file + ! ibufr is the bufr id to be used in subsequent calls call codes_bufr_new_from_file(ifile,ibufr,iret) - do while (iret/=CODES_END_OF_FILE) - write(*,*) 'message: ',count - - ! we need to instruct ecCodes to expand all the descriptors - ! i.e. unpack the data values - call codes_set(ibufr,"unpack",1); - - ! In what follows we rely on the fact that for - ! temperature significant levels the value of key - ! verticalSoundingSignificance is 4 (see flag table 8001 for details). - ! - ! In our BUFR message verticalSoundingSignificance is always followed by - ! geopotential, airTemperature, dewpointTemperature, - ! windDirection, windSpeed and pressure. - ! - ! So in order to access any of these keys we need to use the - ! condition: verticalSoundingSignificance=4. - - ! ---- Get pressure --------------------------- - call codes_get(ibufr,'/verticalSoundingSignificance=4/pressure',presVal); - - ! ---- Get gepotential ------------------------ - call codes_get(ibufr,'/verticalSoundingSignificance=4/geopotential',geoVal) - - ! ---- Get temperature -------------------------------- - call codes_get(ibufr,'/verticalSoundingSignificance=4/airTemperature',tVal) - - ! ---- Get dew point temperature ----------------------- - call codes_get(ibufr,'/verticalSoundingSignificance=4/dewpointTemperature',tdVal) - - ! ---- Check that all arrays are same size - if (size(presVal)/=size(geoVal) .or. size(tVal)/=size(tdVal) .or. size(tdVal)/=size(presVal)) then - print *,'inconsistent array dimension' - exit - endif - num=size(presVal) - - ! ---- Print the values -------------------------------- - write(*,*) 'level pres geo t td' - write(*,*) "--------------------------------------" - - do i=1,num - write(*,*) i,presVal(i),geoVal(i),tVal(i),tdVal(i) - end do - - - ! free arrays - deallocate(presVal) - deallocate(geoVal) - deallocate(tVal) - deallocate(tdVal) - + call codes_set(ibufr,'unpack',1) + call codes_get(ibufr,'timePeriod',timePeriod) + call codes_get(ibufr,'pressure',pressure) + call codes_get(ibufr,'extendedVerticalSoundingSignificance',extendedVerticalSoundingSignificance) + call codes_get(ibufr,'geopotentialHeight',geopotentialHeight) + call codes_get(ibufr,'latitudeDisplacement',latitudeDisplacement) + call codes_get(ibufr,'longitudeDisplacement',longitudeDisplacement) + call codes_get(ibufr,'airTemperature',airTemperature) + call codes_get(ibufr,'dewpointTemperature',dewpointTemperature) + call codes_get(ibufr,'windDirection',windDirection) + call codes_get(ibufr,'windSpeed',windSpeed) + call codes_get(ibufr,'blockNumber',blockNumber) + call codes_get(ibufr,'stationNumber',stationNumber) + print *,'station',blockNumber,stationNumber + print *,'timePeriod pressure geopotentialHeight latitudeDisplacement & + &longitudeDisplacement airTemperature windDirection windSpeed significance' + do i=1,size(windSpeed) + write(*,'(I5,6X,F9.1,2X,F9.2,10X,F8.2,14X,F8.2,16X,F8.2,6X,F8.2,4X,F8.2,4X,I0)') timePeriod(i),pressure(i),& + &geopotentialHeight(i),latitudeDisplacement(i),& + &longitudeDisplacement(i),airTemperature(i),windDirection(i),windSpeed(i),extendedVerticalSoundingSignificance(i) + enddo + ! free arrays + deallocate(timePeriod) + deallocate(pressure) + deallocate(geopotentialHeight) + deallocate(latitudeDisplacement) + deallocate(longitudeDisplacement) + deallocate(airTemperature) + deallocate(dewpointTemperature) + deallocate(windDirection) + deallocate(windSpeed) + deallocate(extendedVerticalSoundingSignificance) ! release the bufr message call codes_release(ibufr) - ! load the next bufr message call codes_bufr_new_from_file(ifile,ibufr,iret) - count=count+1 - - end do - -! close file + end do + ! close file call codes_close_file(ifile) - end program bufr_read_temp diff --git a/examples/F90/bufr_read_temp.sh b/examples/F90/bufr_read_temp.sh index ec84eccf7..737a088c8 100755 --- a/examples/F90/bufr_read_temp.sh +++ b/examples/F90/bufr_read_temp.sh @@ -18,18 +18,15 @@ label="bufr_read_temp_f" fTmp=${label}.tmp.txt rm -f $fTmp | true -#We check "temp_101.bufr". The path is -#hardcoded in the example +# The path to the BUFR file is hardcoded in the example REDIRECT=/dev/null -#Write the key values into a file +# Run the example ${examples_dir}/eccodes_f_bufr_read_temp 2> $REDIRECT > $fTmp #TODO: check the results -#cat $fTmp - #Clean up -rm -f $fTmp | true +rm -f $fTmp diff --git a/examples/python/bufr_read_temp.py b/examples/python/bufr_read_temp.py index cc45e1360..1579290ee 100644 --- a/examples/python/bufr_read_temp.py +++ b/examples/python/bufr_read_temp.py @@ -1,3 +1,4 @@ +# # Copyright 2005-2015 ECMWF. # # This software is licensed under the terms of the Apache Licence Version 2.0 @@ -6,12 +7,11 @@ # In applying this licence, ECMWF does not waive the privileges and immunities # granted to it by virtue of its status as an intergovernmental organisation # nor does it submit to any jurisdiction. - # # Python implementation: bufr_read_temp # -# Description: how to read temperature significant levels from TEMP BUFR -# messages. +# Description: how to read levels from TEMP BUFR messages. +# # # # Please note that TEMP reports can be encoded in various ways in BUFR. @@ -19,94 +19,50 @@ # messages than the one used in the example. It is advised to use bufr_dump to # understand the structure of the messages. # - import traceback import sys - from eccodes import * -INPUT = '../../data/bufr/temp_101.bufr' +INPUT = '../../data/bufr/PraticaTemp.bufr' VERBOSE = 1 # verbose error reporting - def example(): - # open bufr file f = open(INPUT) - cnt = 0 - - # loop for the messages in the file + # loop over the messages in the file while 1: # get handle for message - gid = codes_bufr_new_from_file(f) - if gid is None: + bufr = codes_bufr_new_from_file(f) + if bufr is None: break - print "message: %s" % cnt - # we need to instruct ecCodes to expand all the descriptors - # i.e. unpack the data values - codes_set(gid, 'unpack', 1) - - # In what follows we rely on the fact that for - # temperature significant levels the value of key - # verticalSoundingSignificance is 4 (see flag table 8001 for details). - - # We also make use of the fact that in our BUFR message - # verticalSoundingSignificance is always followed by geopotential, - # airTemperature, dewpointTemperature, - # windDirection, windSpeed and pressure. - - # Get the number of the temperature significant levels. - - # We find out the number of temperature significant levels by - # counting how many pressure values we have on these levels. - - numSigT = codes_get_size( - gid, "/verticalSoundingSignificance=4/pressure") - print ' Number of temperature significant levels %ld' % (numSigT) - - # Get pressure - sigt_pres = codes_get_array( - gid, "/verticalSoundingSignificance=4/pressure") - - # Get gepotential - sigt_geo = codes_get_array( - gid, "/verticalSoundingSignificance=4/geopotential") - - # Get temperature - sigt_t = codes_get_array( - gid, "/verticalSoundingSignificance=4/airTemperature") - - # Get dew point - sigt_td = codes_get_array( - gid, "/verticalSoundingSignificance=4/dewpointTemperature") - - # Check that all arrays are same size - if len(sigt_pres) != numSigT or len(sigt_geo) != numSigT or \ - len(sigt_t) != numSigT or len(sigt_td) != numSigT: - print 'inconsistent array dimension' - return 1 - - # Print the values - print "lev pres geo t td" - print "-------------------------------" - - for i in xrange(numSigT): - print "%3d %6.0f %6.0f %.1f %.1f" % (i + 1, sigt_pres[i], - sigt_geo[i], sigt_t[i], - sigt_td[i]) - + # i.e. unpack the data section + codes_set(bufr, 'unpack', 1) + # get all the timePeriods + timePeriod = codes_get_array(bufr, "timePeriod") + pressure = codes_get_array(bufr, "pressure") + extendedVerticalSoundingSignificance = codes_get_array(bufr, "extendedVerticalSoundingSignificance") + geopotentialHeight = codes_get_array(bufr, "geopotentialHeight") + latitudeDisplacement = codes_get_array(bufr, "latitudeDisplacement") + longitudeDisplacement = codes_get_array(bufr, "longitudeDisplacement") + airTemperature = codes_get_array(bufr, "airTemperature") + dewpointTemperature = codes_get_array(bufr, "dewpointTemperature") + windDirection = codes_get_array(bufr, "windDirection") + windSpeed = codes_get_array(bufr, "windSpeed") + blockNumber = codes_get(bufr, "blockNumber") + stationNumber = codes_get(bufr, "stationNumber") + print 'station %d%d' % (blockNumber,stationNumber) + print 'timePeriod pressure geopotentialHeight latitudeDisplacement longitudeDisplacement airTemperature windDirection windSpeed significance' + for i in range(0,len(windSpeed)-1): + print timePeriod[i],pressure[i],geopotentialHeight[i],latitudeDisplacement[i],longitudeDisplacement[i],airTemperature[i],windDirection[i],windSpeed[i],extendedVerticalSoundingSignificance[i] cnt += 1 - # delete handle - codes_release(gid) - + codes_release(bufr) # close the file f.close() - def main(): try: example() @@ -115,8 +71,7 @@ def main(): traceback.print_exc(file=sys.stderr) else: print >>sys.stderr, err.msg - return 1 - if __name__ == "__main__": sys.exit(main()) + diff --git a/examples/python/bufr_read_temp.sh b/examples/python/bufr_read_temp.sh index 794b64efe..5af67aaba 100755 --- a/examples/python/bufr_read_temp.sh +++ b/examples/python/bufr_read_temp.sh @@ -16,20 +16,14 @@ label="bufr_read_temp_p" #Define tmp file fTmp=${label}.tmp.txt -rm -f $fTmp | true - -#We check "temp_101.bufr". The path is -#hardcoded in the example +rm -f $fTmp REDIRECT=/dev/null -#Write the key values into a file +# Run it. The path to the BUFR file is hardcoded in the example $PYTHON $examples_src/bufr_read_temp.py 2> $REDIRECT > $fTmp #TODO: check the results -#cat $fTmp - #Clean up rm -f $fTmp -