Bug #64022 for Geo-ReadGRIB: Latitude Increment (LiInc) problem with German's Weather Service (DWD) GRIBs

Fri Dec 17 21:39:20 2010 frank.l.cox [...] gmail.com - Ticket created

Subject:	Latitude Increment (LiInc) problem with German's Weather Service (DWD) GRIBs
From:	frank.l.cox [...] gmail.com

Jose Manuel Fernandez Jaramillo emailed: Hi Frank, First thank you for your perl module for reading grib files, It is helping me to clean some dirty scripts to get the data. My problem is that I am trying to read some files from the German's Weather Service (DWD) and when they write the resolution they have 0.063, and It should be 0.0625. Your PERL module uses the 0.063 and then as I have a 241x417 grid, there is an important shift in the coordinates. According with your script I am able to get only a 239x413 grid. And in a first hacking I found that the temporal file that your script creates has all the elements that are supposed in order to fill the complete grid. Reading the documentation I found a solution forcing the increment by $w->LaInc(0.0625), But as the same time you say not to use in production code. Is there any alternative that you have or planning to have that I can use? may be you can read the number of cells from wgrib and compare the consistency of the increment value. thank you for your attention, Jose Manuel P.D. this is the verbose output for one record in my Grib file: rec 226:48076650:date 2010121700 T kpds5=11 kpds6=105 kpds7=2 levels=(0,2) grid=255 2 m above gnd 78hr fcst: bitmap: 0 undef T=temperature [K] timerange 0 P1 78 P2 0 TimeU 1 nx 417 ny 241 GDS grid 0 num_in_ave 0 missing 0 center 78 subcenter 255 process 135 Table 2 scan: WE:NS winds(N/S) latlon: lat 62.500000 to 47.500000 by 0.063000 nxny 100497 long -12.500000 to 13.500000 by 0.063000, (417 x 241) scan 0 mode 128 bdsgrid 1 min/max data 253.937 283.613 num bits 16 BDS_Ref 253.937 DecScale 0 BinScale -11

Fri Dec 17 21:50:03 2010 frank.l.cox [...] gmail.com - Cc josemanuelfernandezjaramillo@yahoo.com added

Fri Dec 17 22:00:13 2010 frank.l.cox [...] gmail.com - Correspondence added

From the description it looks like a problem with DWD GRIBs or the way Geo::ReadGRIB handles them. I'll get some more details from the reporter and get some of these files and try to reproduce the problem... On Fri Dec 17 21:39:20 2010, FRANKCOX wrote: Show quoted text

> Jose Manuel Fernandez Jaramillo emailed: > > Hi Frank, > > First thank you for your perl module for reading grib files, It is > helping me to > clean some dirty scripts to get the data. > > My problem is that I am trying to read some files from the German's > Weather > Service (DWD) and when they write the resolution they have 0.063, and > It should > be 0.0625. Your PERL module uses the 0.063 and then as I have a > 241x417 grid, > there is an important shift in the coordinates. According with your > script I am > able to get only a 239x413 grid. And in a first hacking I found that > the > temporal file that your script creates has all the elements that are > supposed in > order to fill the complete grid. > > Reading the documentation I found a solution forcing the increment by > $w->LaInc(0.0625), But as the same time you say not to use in > production code. > > Is there any alternative that you have or planning to have that I can > use? > > may be you can read the number of cells from wgrib and compare the > consistency > of the increment value. > > thank you for your attention, Jose Manuel > > P.D. > > this is the verbose output for one record in my Grib file: > > > rec 226:48076650:date 2010121700 T kpds5=11 kpds6=105 kpds7=2 > levels=(0,2) > grid=255 2 m above gnd 78hr fcst: bitmap: 0 undef > T=temperature [K] > timerange 0 P1 78 P2 0 TimeU 1 nx 417 ny 241 GDS grid 0 num_in_ave > 0 missing > 0 > center 78 subcenter 255 process 135 Table 2 scan: WE:NS winds(N/S) > latlon: lat 62.500000 to 47.500000 by 0.063000 nxny 100497 > long -12.500000 to 13.500000 by 0.063000, (417 x 241) scan 0 > mode 128 > bdsgrid 1 > min/max data 253.937 283.613 num bits 16 BDS_Ref 253.937 DecScale > 0 > BinScale -11

Fri Dec 17 22:00:13 2010 frank.l.cox [...] gmail.com - Status changed from 'new' to 'open'

Fri Dec 17 22:01:22 2010 frank.l.cox [...] gmail.com - Severity Normal added

Fri Dec 17 22:01:22 2010 frank.l.cox [...] gmail.com - Broken in 1.4 added

Mon Dec 20 13:16:39 2010 frank.l.cox [...] gmail.com - Correspondence added

José Manuel sent a Perl script and a sample DWD GRIB file that illustrates the problem he reported. I'm attaching the script. The GRIB file is too big to upload. I'll be looking at this problem as soon as I can. --Frank Cox

Subject:

readgrib_dwd_test.pl

#!/usr/bin/perl -w use Geo::ReadGRIB; $timeindex = 0; open FILE1, ">lon.txt" or die $!; open FILE2, ">lat.txt" or die $!; open FILE3, ">val.txt" or die $!; $w = new Geo::ReadGRIB "dwd_lm2_2007120100.gb1"; $w->getFullCatalog(); print $w->show(); #print $w->getParam(show)."\n"; $wtime = $w->getParam(TIME); $wlasttime = $w->getParam(LAST_TIME); $wla1 = $w->getParam(La1); $wla2 = $w->getParam(La2); $wlainc = $w->getParam(LaInc); $wlo1 = $w->getParam(Lo1); $wlo2 = $w->getParam(Lo2); $wloinc = $w->getParam(LoInc); $wfilename = $w->getParam(fileName); print $wtime," ",$wlasttime," ",$wla1," ",$wla2," ",$wlainc," ",$wlo1," ",$wlo2," ",$wloinc," ",$wfilename; #$w->LaInc(0.0625); #$w->LoInc(0.0625); print "\nlainc new: ".$w->getParam(LaInc)."\n"; for ( $timedata = $wtime + $timeindex * 60 * 60 ; $timedata <= $wtime + $timeindex * 60 * 60 ; #$timedata <= $wlasttime ; $timedata = $timedata + 60 * 60 ) { $split = $w->extractLaLo( [ 'PRMSL', 'TMP', 'UGRD', 'VGRD' ], $wla1, $wlo1, $wla2, $wlo2, $timedata ); #Geo::ReadGRIB::PlaceIterator #die $w->getError if $w->getError; #print "$w->getError\n"; if ( $w->getError ) { $split = $w->extractLaLo( ['TMP'], $wla1, $wlo1, $wla2, $wlo2, $timedata ) ; #Geo::ReadGRIB::PlaceIterator } else { print "ok\n"; } if ( $w->getError ) { die $w->getError } else { print "ok\n"; } $vindex = 0; my @vtmp; my @vugrd; my @vvgrd; my @vprmsl; print "numlat: ".$split->numLat." numlong: ".$split->numLong."\n"; for my $y ( 0 .. $split->numLat - 1 ) { for my $x ( 0 .. $split->numLong - 1 ) { my $place = $split->current; # $place is a Geo::ReadGRIB::Place object $time = $place->thisTime; $latitude = $place->lat; $longitude = $place->long; $data_types = $place->types; # an array ref of type names #print "yx $y $x $latitude $longitude\n"; @data_deferenced = @$data_types; #print "numero elementos ".($#data_deferenced)."\n"; #print "@data_deferenced,\n"; #print ("($latitude $longitude) "); for ( $param = 0 ; $param <= $#data_deferenced ; $param++ ) { $data = $place->data( $$data_types[$param] ); if ( $$data_types[$param] eq "TMP" ) { $vtmp[$vindex] = $data; print FILE1 "$longitude "; print FILE2 "$latitude "; print FILE3 "$data "; } elsif ( $$data_types[$param] eq "UGRD" ) { $vugrd[$vindex] = $data; } elsif ( $$data_types[$param] eq "VGRD" ) { $vvgrd[$vindex] = $data; } elsif ( $$data_types[$param] eq "PRMSL" ) { $vprmsl[$vindex] = $data; } else { print "Data not stored\n"; } } $split->next; $vindex++; } print FILE1 "\n"; print FILE2 "\n"; print FILE3 "\n"; } print "vindex: $vindex, Last indexvector result: $#vtmp\n"; } exit 0;

Sun May 19 20:55:11 2013 frank.l.cox [...] gmail.com - Correspondence added

grib2 files are still not supported

Sun May 19 20:55:11 2013 frank.l.cox [...] gmail.com - Status changed from 'open' to 'resolved'