diff --git a/docs/web_docs/set_grib_max_bits.html b/docs/web_docs/set_grib_max_bits.html index d3bdd36f..c81270c7 100644 --- a/docs/web_docs/set_grib_max_bits.html +++ b/docs/web_docs/set_grib_max_bits.html @@ -133,6 +133,30 @@

Introduction

to 25 after accounting for the limits of the packing and unpacking routines. +

The complex-2 and complex-3 compression schemes also +store scaled integers but different. The complex-1 stores +the scaled intger values, so it is limited to a 25 bit +integers. Hewever, complex-2 stores the delta (the difference +from previously defined grid point). The magnitude of delta +is limited by the original series but the delta can be +either positive or negative. So under worst case, the +deltas need one more bit to store the scaled integers. +So if you scale the grid point values to 24 bits, you +will have a scaled deltas being 25 bits. The complex-3 +packing stores the deltas of the delta series. This is +good for smooth fields, but for the worst case could +cause two extra bits to store the original grid data. +So at worst case, the data would need to be scaled at 23 +bits. + +

The extra bits needed for complex-2 and complex-3 are +very rare. It was first encountered in 11/2024 which +is about a decade after the complex packing was introduced. +In 11/2024, code was added to the complex packing to detect +the increase in bits needed and reduce the precision if it +cause the precision to be greater than 25 bits. + +

Grib fields are usually stored as a scaled integers, and usually you don't need 25 bits of precision. For example, temperature to the nearest 0.1 degree, means diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt index 37b853fc..4bcd1879 100644 --- a/tests/CMakeLists.txt +++ b/tests/CMakeLists.txt @@ -109,6 +109,7 @@ copy_test_data(ref_simple_packing.grib2.spread.txt) copy_test_data(ref_new_grid_gdt_32769.grib2) copy_test_data(png_4bits.png) copy_test_data(large_png.grb2) +copy_test_data(ref_c3_overflow.grib2) # Run these shell tests. shell_test(run_wgrib2_tests) diff --git a/tests/data/ref_c3_overflow.grib2 b/tests/data/ref_c3_overflow.grib2 new file mode 100644 index 00000000..612c3cb7 Binary files /dev/null and b/tests/data/ref_c3_overflow.grib2 differ diff --git a/tests/run_wgrib2_tests.sh b/tests/run_wgrib2_tests.sh index e1252285..127cdac0 100755 --- a/tests/run_wgrib2_tests.sh +++ b/tests/run_wgrib2_tests.sh @@ -145,7 +145,7 @@ if [ `echo "$line" | grep -c ":rpn_rms=0:"` -ne 1 ] ; then exit 1 fi - +../wgrib2/wgrib2 data/ref_c3_overflow.grib2 -set_grib_type c3 -set_grib_max_bits 25 -grib_out c3_overflow.grib2 echo "*** SUCCESS!" exit 0 diff --git a/wgrib2/complex_pk.c b/wgrib2/complex_pk.c index dc3e5fd1..e1c3af72 100644 --- a/wgrib2/complex_pk.c +++ b/wgrib2/complex_pk.c @@ -84,13 +84,13 @@ static int sizeofsection2(int mn, int mx, int n, int ref_bits, int width_bits, i } static int size_all(struct section *s, int ref_bits, int width_bits, int has_undef) { - unsigned int bytes, bits; + unsigned int bytes, bits; bytes = bits = 0; while (s) { bits += sizeofsection(s, ref_bits, width_bits, has_undef); - bytes += (bits / 8); - bits = bits % 8; + bytes += (bits / 8); + bits = bits % 8; s = s->tail; } return (int) (bytes + (bits != 0)); @@ -722,8 +722,7 @@ int complex_grib_out(unsigned char **sec, float *data, unsigned int ndata, if (int_min_max_array(v, nndata, &vmn, &vmx) != 0) fatal_error("complex_pk: int_min_max error",""); } -//fprintf(stderr , "\n pre process v[i] data extri_0 %d extra_1 %d\n",extra_0, extra_1); -//for (i = 0; i < nndata;i++) { +//fprintf(stderr , "\n pre process v[i] data extri_0 %d extra_1 %d\n",extra_0, extra_1); //for (i = 0; i < nndata;i++) { //fprintf(stderr," %d:%d ", i, v[i]); //} //fprintf(stderr,"\n"); @@ -732,6 +731,25 @@ int complex_grib_out(unsigned char **sec, float *data, unsigned int ndata, printf("2: vmx %d vmn %d nbits %d\n", vmx, vmn, find_nbits(vmx-vmn+has_undef)); #endif + /* + c2: values are deltas, which are both plus and minus + so range of values is up to twice range of input (one more bit) + c3: values are deltas of the deltas + the range of deltas are up to 2x initial values + the potential rannge of deltas of deltas are 4x initial values + */ + + j = find_nbits(vmx-vmn+has_undef); + if (j > 25) { + free(v); + free(sec5); + free(sec6); + max_bits -= 1; + if (packing_mode == 3) max_bits -= 1; + return complex_grib_out(sec, data, ndata, use_scale, dec_scale, bin_scale, + wanted_bits, max_bits, packing_mode, use_bitmap, out); + } + if (has_undef) { #pragma omp parallel for private(i) for (i = 0; i < nndata; i++) {