Re: [netcdf-java] Grib1 Gaussian Grid Relative Humidity Interpolation Issue

Hi,

hmm… did you have the chance to look deeper into this issue, yet?

Thanks,
Jochen

Am 25.11.2013 um 17:35 schrieb Jochen Kähler <jkaehler@xxxxxxxxxxxxxxx>:

> Hi John,
> 
> that an easy one: let’s go for option #4! ;-) ...
> 
> Well, regarding to our current project, linear interpolation would be totally 
> fine. If we could get some configuration switch or even some option to change 
> interpolationmethod programmatically, that would be excellent.
> 
> Wgrib provides access to the original numbers - that would be an interesting 
> option as well. As my colleague Martijn already pointed out here:  
> http://www.unidata.ucar.edu/mailing_lists/archives/netcdf-java/2013/msg00125.html
>  , we also have some issues around the chosen resolution, which we could work 
> around in a somewhat more elegant way if we had access to the raw data...
> 
> Jochen
> 
> Am 19.11.2013 um 18:37 schrieb John Caron <caron@xxxxxxxxxxxxxxxx>:
> 
>> Hi Jochen:
>> 
>> We call these "reduced grids" or "thin grids", their use is technically 
>> independent of whether its a gaussian lat/lon (just means the lats are 
>> spaced with a gaussian), although in practice they often go together 
>> (perhaps always in the model - im just talking about the output files).  We 
>> probably stole the cubic spline code from GEMPAK/nawips which made us think 
>> it was a reassonable thing to do.
>> 
>> IFAIU, cubic splines allow the interpolated numbers to be outside the range 
>> of existing ones, so at this point i dont suspect an error.
>> 
>> You are right that we could use the actual min, max of the numbers to clip 
>> the cubic spline results.
>> 
>> Not sure what the right thing to do here is. we could clip it based on 
>> actual min/max (easy). We could try to make an option to allow linear 
>> interpolation, which would prevent the interpolated numbers from exceeding 
>> the min/max of existing numbers (harder). We could try to allow user access 
>> to the original numbers (harder, need lower level API). We could construct a 
>> quantity X conserving algorithm, which reads the users' mind as to what X 
>> should be (requires faster than light technology).
>> 
>> What does wgrib do?
>> 
>> Suggestions from anyone are welcome.
>> 
>> John
>> 
>>>> 
>>>> On 11/15/2013 6:04 AM, Jochen Kähler wrote:
>>>>> Hi all,
>>>>> 
>>>>> we’re implementing some service that reads ECMWF 0.125° Grib files in
>>>>> gaussian grid. The NetCDF-Library automatically interpolates the
>>>>> gaussian grid to a regular lat/lon-grid. For the relative humidity we
>>>>> found some strange values in output after that interpolation.
>>>>> 
>>>>> Using wgrib I’ve extracted some isobaric relative humidity record to a
>>>>> single grib file called „dpd_test.grib“ (which i could send as well if
>>>>> needed). According to wgrib, it is an gaussian grid and min and max
>>>>> values are 0% and 100% - as I would expect them for a relative humidity
>>>>> (see attached output).
>>>>> 
>>>>> Next I’ve written a small program which extracts that relative humidity
>>>>> parameter via netcdf-java-lib and step thru all values to find minimum
>>>>> and maximum in it. Source is also attached - just in case I did
>>>>> something wrong there.
>>>>> 
>>>>> The output of that states the following...
>>>>> 2013-11-14 13:45:15,378 INFO : App -
>>>>> dpd_test.grib/Relative_humidity_isobaric/850.0hPa:
>>>>> min=-4.462653160095215 max=108.53846740722656
>>>>> …indicating the minimum relative humidity is around -4.46% and the
>>>>> maximum at 108.54%.
>>>>> 
>>>>> Browsing your code on Github I think you’re using cubic interpolation to
>>>>> convert gaussian to regular grid, which might cause this problem.
>>>>> https://github.com/Unidata/thredds/blob/6b1052455e597797f3a5980165a23292172920ce/grib/src/main/java/ucar/nc2/grib/QuasiRegular.java
>>>>>  Line
>>>>> 219 - 221:
>>>>> outpt[oIdx] = (float) (a * inpt[iIdx + low] + b * inpt[iIdx + hi]
>>>>>         + ((a * a * a - a) * y2d[low]
>>>>>         + (b * b * b - b) * y2d[hi]) / 6.0);
>>>>> Is the problem in the parametrization of this formula?
>>>>> 
>>>>> Is it possible to use Linear Interpolation here?
>>>>> 
>>>>> Kind Regards,
>>>>> Jochen
>>>>> 
>>>>> 
>>>>> 
>>>>> — wgrib output ---------------
>>>>> Losty-MacBook:Desktop jkaehler$ ~/wgrib/wgrib -V dpd_test.grib
>>>>> rec 1:0:date 2013111300 R kpds5=157 kpds6=100 kpds7=850 levels=(3,82)
>>>>> grid=255 850 mb 6hr fcst:
>>>>>   R=Relative humidity [%]
>>>>>   timerange 0 P1 6 P2 0 TimeU 1  nx -1 ny 800 GDS grid 4 num_in_ave 0
>>>>> missing 0
>>>>>   center 98 subcenter 0 process 143 Table 128 scan: WE:NS winds(N/S)
>>>>>   thinned gaussian: lat  89.828000 to -89.828000
>>>>>           long 0.000000 to 359.900000, 843490 grid pts   (-1 x 800)
>>>>> scan 0 mode 0 bdsgrid 1
>>>>>       18   25   32   40   45   50   60   60   72   72   75   81   90
>>>>> 96  100
>>>>>      108  120  120  125  128  144  144  150  160  160  180  180  192
>>>>>  192  200
>>>>>      200  216  216  225  240  240  240  250  250  256  270  288  288
>>>>>  288  300
>>>>>      300  320  320  320  324  360  360  360  360  360  360  375  375
>>>>>  384  400
>>>>>      400  400  405  432  432  432  432  450  450  450  480  480  480
>>>>>  480  480
>>>>>      486  500  500  512  512  540  540  540  540  540  576  576  576
>>>>>  576  576
>>>>>      576  600  600  600  600  640  640  640  640  640  640  640  648
>>>>>  675  675
>>>>>      675  675  675  720  720  720  720  720  720  720  729  729  750
>>>>>  750  750
>>>>>      750  768  768  768  800  800  800  800  800  800  810  864  864
>>>>>  864  864
>>>>>      864  864  864  864  864  864  900  900  900  900  900  900  900
>>>>>  960  960
>>>>>      960  960  960  960  960  960  960  960  960  960  972  972 1000
>>>>> 1000 1000
>>>>>     1000 1000 1000 1024 1024 1024 1024 1024 1080 1080 1080 1080 1080
>>>>> 1080 1080
>>>>>     1080 1080 1080 1080 1125 1125 1125 1125 1125 1125 1125 1125 1125
>>>>> 1152 1152
>>>>>     1152 1152 1152 1152 1200 1200 1200 1200 1200 1200 1200 1200 1200
>>>>> 1200 1200
>>>>>     1215 1215 1215 1215 1280 1280 1280 1280 1280 1280 1280 1280 1280
>>>>> 1280 1280
>>>>>     1280 1280 1280 1280 1280 1296 1296 1296 1296 1350 1350 1350 1350
>>>>> 1350 1350
>>>>>     1350 1350 1350 1350 1350 1350 1350 1350 1350 1440 1440 1440 1440
>>>>> 1440 1440
>>>>>     1440 1440 1440 1440 1440 1440 1440 1440 1440 1440 1440 1440 1440
>>>>> 1440 1440
>>>>>     1440 1440 1440 1440 1440 1440 1440 1458 1458 1458 1458 1458 1458
>>>>> 1458 1500
>>>>>     1500 1500 1500 1500 1500 1500 1500 1500 1500 1500 1500 1500 1500
>>>>> 1500 1500
>>>>>     1500 1536 1536 1536 1536 1536 1536 1536 1536 1536 1536 1536 1536
>>>>> 1536 1536
>>>>>     1536 1536 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600
>>>>> 1600 1600
>>>>>     1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600
>>>>> 1600 1600
>>>>>     1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600
>>>>> 1600 1600
>>>>>     1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600
>>>>> 1600 1600
>>>>>     1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600
>>>>> 1600 1600
>>>>>     1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600
>>>>> 1600 1600
>>>>>     1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600
>>>>> 1600 1600
>>>>>     1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600
>>>>> 1600 1600
>>>>>     1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600
>>>>> 1600 1600
>>>>>     1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600
>>>>> 1600 1600
>>>>>     1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600 1600
>>>>> 1600 1600
>>>>>     1600 1600 1600 1536 1536 1536 1536 1536 1536 1536 1536 1536 1536
>>>>> 1536 1536
>>>>>     1536 1536 1536 1536 1500 1500 1500 1500 1500 1500 1500 1500 1500
>>>>> 1500 1500
>>>>>     1500 1500 1500 1500 1500 1500 1458 1458 1458 1458 1458 1458 1458
>>>>> 1440 1440
>>>>>     1440 1440 1440 1440 1440 1440 1440 1440 1440 1440 1440 1440 1440
>>>>> 1440 1440
>>>>>     1440 1440 1440 1440 1440 1440 1440 1440 1440 1440 1440 1350 1350
>>>>> 1350 1350
>>>>>     1350 1350 1350 1350 1350 1350 1350 1350 1350 1350 1350 1296 1296
>>>>> 1296 1296
>>>>>     1280 1280 1280 1280 1280 1280 1280 1280 1280 1280 1280 1280 1280
>>>>> 1280 1280
>>>>>     1280 1215 1215 1215 1215 1200 1200 1200 1200 1200 1200 1200 1200
>>>>> 1200 1200
>>>>>     1200 1152 1152 1152 1152 1152 1152 1125 1125 1125 1125 1125 1125
>>>>> 1125 1125
>>>>>     1125 1080 1080 1080 1080 1080 1080 1080 1080 1080 1080 1080 1024
>>>>> 1024 1024
>>>>>     1024 1024 1000 1000 1000 1000 1000 1000  972  972  960  960  960
>>>>>  960  960
>>>>>      960  960  960  960  960  960  960  900  900  900  900  900  900
>>>>>  900  864
>>>>>      864  864  864  864  864  864  864  864  864  810  800  800  800
>>>>>  800  800
>>>>>      800  768  768  768  750  750  750  750  729  729  720  720  720
>>>>>  720  720
>>>>>      720  720  675  675  675  675  675  648  640  640  640  640  640
>>>>>  640  640
>>>>>      600  600  600  600  576  576  576  576  576  576  540  540  540
>>>>>  540  540
>>>>>      512  512  500  500  486  480  480  480  480  480  450  450  450
>>>>>  432  432
>>>>>      432  432  405  400  400  400  384  375  375  360  360  360  360
>>>>>  360  360
>>>>>      324  320  320  320  300  300  288  288  288  270  256  250  250
>>>>>  240  240
>>>>>      240  225  216  216  200  200  192  192  180  180  160  160  150
>>>>>  144  144
>>>>>      128  125  120  120  108  100   96   90   81   75   72   72   60
>>>>> 60   50
>>>>>       45   40   32   25   18
>>>>>   min/max data 0 100  num bits 8  BDS_Ref 0  DecScale 0 BinScale -1
>>>>> ——
>>>>> 
>>>>> 
>>>>> — Java-Source ——
>>>>> GridDataset dataset = GridDataset.open(gribFile.getAbsolutePath());
>>>>> 
>>>>> try {
>>>>> GridDatatype grid = 
>>>>> dataset.findGridDatatype("Relative_humidity_isobaric");
>>>>> if(grid != null) {
>>>>> CoordinateAxis1D isbl = grid.getCoordinateSystem().getVerticalAxis();
>>>>> double[] isblValues = isbl.getCoordValues();
>>>>> for (int isblIdx = 0; isblIdx < isblValues.length; isblIdx++) {
>>>>> double[] data = (double[]) grid.readDataSlice(0, isblIdx, -1,
>>>>> -1).get1DJavaArray(double.class);
>>>>> double min = Double.MAX_VALUE;
>>>>> double max = Double.MIN_VALUE;
>>>>> for (double d : data) {
>>>>> min = Math.min(min, d);
>>>>> max = Math.max(max, d);
>>>>>        }
>>>>> LOG.info("{}/{}/{}{}: min={} max={}", gribFile.getName(),
>>>>> grid.getFullName(), isblValues[isblIdx], isbl.getUnitsString(), min, max);
>>>>> }
>>>>> }
>>>>> } finally {
>>>>> if (dataset != null) {
>>>>> dataset.close();
>>>>> }
>>>>> }
>>>>> 
>>>>> --------
>>>>> 
>>>>> 
>>>>> _______________________________________________
>>>>> netcdf-java mailing list
>>>>> netcdf-java@xxxxxxxxxxxxxxxx
>>>>> For list information or to unsubscribe, visit: 
>>>>> http://www.unidata.ucar.edu/mailing_lists/
>>>>> 
>>>> 
>>>> _______________________________________________
>>>> netcdf-java mailing list
>>>> netcdf-java@xxxxxxxxxxxxxxxx
>>>> For list information or to unsubscribe, visit: 
>>>> http://www.unidata.ucar.edu/mailing_lists/
>>> 
>> 
>> _______________________________________________
>> netcdf-java mailing list
>> netcdf-java@xxxxxxxxxxxxxxxx
>> For list information or to unsubscribe, visit: 
>> http://www.unidata.ucar.edu/mailing_lists/ 
> 

Attachment: smime.p7s
Description: S/MIME cryptographic signature

  • 2013 messages navigation, sorted by:
    1. Thread
    2. Subject
    3. Author
    4. Date
    5. ↑ Table Of Contents
  • Search the netcdf-java archives: