diff --git a/src/grib_accessor_class_from_scale_factor_scaled_value.c b/src/grib_accessor_class_from_scale_factor_scaled_value.c index bf339246d..70d7d9763 100644 --- a/src/grib_accessor_class_from_scale_factor_scaled_value.c +++ b/src/grib_accessor_class_from_scale_factor_scaled_value.c @@ -180,6 +180,79 @@ static double eval_value_factor(int64_t value, int64_t factor) return (double)value * pow(10.0, -factor); } +/* Return 0 if no error and set the two outputs 'value' and 'factor' + * value cannot exceed maximum_value and factor cannot exceed maximum_factor + */ +static int get_scaled_value_and_scale_factor_algorithm1( + double input, int64_t maximum_value, int64_t maximum_factor, + int64_t* ret_value, int64_t* ret_factor) +{ + int64_t factor = 0; + int64_t value = 0; + + factor = floor(log10(maximum_value)) - floor(log10(input < 0 ? -input : input)); + if (factor > maximum_factor) { + return GRIB_INTERNAL_ERROR; + } + value = (int64_t)round(input * pow(10, factor)); + while ((value % 10 == 0) && (factor > 0)) { + value /= 10; + factor--; + } + if (value >= maximum_value) + return GRIB_INTERNAL_ERROR; + if (factor >= maximum_factor) + return GRIB_INTERNAL_ERROR; + + *ret_factor = factor; + *ret_value = value; + return GRIB_SUCCESS; +} + +/* Return 0 if no error and set the two outputs 'value' and 'factor' + * value cannot exceed maximum_value and factor cannot exceed maximum_factor + */ +static int get_scaled_value_and_scale_factor_algorithm2( + double input, int64_t maximum_value, int64_t maximum_factor, + int64_t* ret_value, int64_t* ret_factor) +{ + int64_t factor = 0, prev_factor = 0; + int64_t value = 0, prev_value = 0; + double exact = input; + const float epsilon = float_epsilon(); + int is_negative = 0; + + /* Loop until we find a close enough approximation. Keep the last good values */ + if (exact < 0) { + is_negative = 1; + exact *= -1; + } + factor = prev_factor = 0; + value = prev_value = round(exact); + while (!is_approximately_equal(exact, eval_value_factor(value, factor), epsilon) && + value < maximum_value && + factor < maximum_factor) { + value = round(exact * pow(10., ++factor)); + if (value > maximum_value || factor > maximum_factor) { + /* One or more maxima exceeded. So stop and use the previous values */ + value = prev_value; + factor = prev_factor; + break; + } + prev_factor = factor; + prev_value = value; + } + + if (is_negative) { + value *= -1; + } + + *ret_factor = factor; + *ret_value = value; + + return GRIB_SUCCESS; +} + static int pack_double(grib_accessor* a, const double* val, size_t* len) { /* See ECC-979 */ @@ -187,12 +260,10 @@ static int pack_double(grib_accessor* a, const double* val, size_t* len) grib_accessor_from_scale_factor_scaled_value* self = (grib_accessor_from_scale_factor_scaled_value*)a; grib_handle* hand = grib_handle_of_accessor(a); int ret = 0; - int64_t factor = 0, prev_factor = 0; - int64_t value = 0, prev_value = 0; + int64_t factor = 0; + int64_t value = 0; double exact = *val; /*the input*/ - const float epsilon = float_epsilon(); - int is_negative = 0; - unsigned long maxval_value, maxval_factor; /*maximum allowable values*/ + int64_t maxval_value, maxval_factor; /*maximum allowable values*/ grib_accessor *accessor_factor, *accessor_value; if (exact == 0) { @@ -220,29 +291,13 @@ static int pack_double(grib_accessor* a, const double* val, size_t* len) maxval_value = (1UL << (accessor_value->length * 8)) - 2; /* exclude missing */ maxval_factor = (1UL << (accessor_factor->length * 8)) - 2; /* exclude missing */ - /* Loop until we find a close enough approximation. Keep the last good values */ - if (exact < 0) { - is_negative = 1; - exact *= -1; - } - factor = prev_factor = 0; - value = prev_value = round(exact); - while (!is_approximately_equal(exact, eval_value_factor(value, factor), epsilon) && - value < maxval_value && - factor < maxval_factor) { - value = round(exact * pow(10., ++factor)); - if (value > maxval_value || factor > maxval_factor) { - /* One or more maxima exceeded. So stop and use the previous values */ - value = prev_value; - factor = prev_factor; - break; + ret = get_scaled_value_and_scale_factor_algorithm1(exact, maxval_value, maxval_factor, &value, &factor); + if (ret) { + ret = get_scaled_value_and_scale_factor_algorithm2(exact, maxval_value, maxval_factor, &value, &factor); + if (ret) { + grib_context_log(a->context, GRIB_LOG_ERROR, "Failed to compute %s and %s from %g", self->scaleFactor, self->scaledValue, exact); + return ret; } - prev_factor = factor; - prev_value = value; - } - - if (is_negative) { - value *= -1; } if ((ret = grib_set_long_internal(hand, self->scaleFactor, factor)) != GRIB_SUCCESS) diff --git a/tests/grib_ecc-1150.sh b/tests/grib_ecc-1150.sh index 83fc00d59..d03dcfede 100755 --- a/tests/grib_ecc-1150.sh +++ b/tests/grib_ecc-1150.sh @@ -88,25 +88,25 @@ ${tools_dir}/grib_set -s lowerLimit=-6.6,upperLimit=-1.02 $tempGrib $temp2 grib_check_key_equals $temp2 scaleFactorOfLowerLimit,scaledValueOfLowerLimit,lowerLimit "1 -66 -6.6" grib_check_key_equals $temp2 scaleFactorOfUpperLimit,scaledValueOfUpperLimit,upperLimit "2 -102 -1.02" -# input factor value -set_lowerLimit_and_check 550 0 550 -set_lowerLimit_and_check -99 0 -99 -set_lowerLimit_and_check 6.77 2 677 -set_lowerLimit_and_check 0.001 3 1 -set_lowerLimit_and_check -6.6 1 -66 -set_lowerLimit_and_check -1.02 2 -102 -set_lowerLimit_and_check 3e-05 5 3 -set_lowerLimit_and_check -3.9e-05 6 -39 - -#set_lowerLimit_and_check 3.14e-06 6 314 -#set_lowerLimit_and_check 3.14e-07 7 314 -#set_lowerLimit_and_check 3e-08 8 3 -#set_lowerLimit_and_check 1.0e-10 10 1 -#set_lowerLimit_and_check 0.03e-06 8 3 -#set_lowerLimit_and_check 3.14e-09 11 314 -#set_lowerLimit_and_check -3.1456e-09 13 -31456 -#set_lowerLimit_and_check 0.0000123 7 123 - +# input factor value +set_lowerLimit_and_check 0 0 0 +set_lowerLimit_and_check 1 0 1 +set_lowerLimit_and_check 550 0 550 +set_lowerLimit_and_check -99 0 -99 +set_lowerLimit_and_check 6.77 2 677 +set_lowerLimit_and_check 0.001 3 1 +set_lowerLimit_and_check -6.6 1 -66 +set_lowerLimit_and_check -1.02 2 -102 +set_lowerLimit_and_check 3e-05 5 3 +set_lowerLimit_and_check -3.9e-05 6 -39 +set_lowerLimit_and_check 3.14e-06 8 314 +set_lowerLimit_and_check 3.14e-07 9 314 +set_lowerLimit_and_check 3e-08 8 3 +set_lowerLimit_and_check -3e-08 8 -3 +set_lowerLimit_and_check 1e-10 10 1 +set_lowerLimit_and_check 3.14e-09 11 314 +set_lowerLimit_and_check 1.23e-05 7 123 +set_lowerLimit_and_check -3.1456e-09 13 -31456 # Clean up rm -f $tempGrib $tempFilt $temp2