diff --git a/VERSION b/VERSION index dc1591835..78d1d9cf4 100644 --- a/VERSION +++ b/VERSION @@ -1 +1 @@ -2.34.0 +2.34.1 diff --git a/definitions/grib2/cfVarName.def b/definitions/grib2/cfVarName.def index 3a2f1f247..889c839e3 100644 --- a/definitions/grib2/cfVarName.def +++ b/definitions/grib2/cfVarName.def @@ -7841,7 +7841,9 @@ parameterCategory = 4 ; parameterNumber = 27 ; typeOfFirstFixedSurface = 168 ; - typeOfSecondFixedSurface = 168 ; + typeOfSecondFixedSurface = 255 ; + scaledValueOfSecondFixedSurface = missing() ; + scaleFactorOfSecondFixedSurface = missing() ; } #Sea water potential temperature tendency due to newtonian relaxation 'thetaodmp' = { @@ -9028,7 +9030,9 @@ parameterCategory = 4 ; parameterNumber = 27 ; typeOfFirstFixedSurface = 168 ; - typeOfSecondFixedSurface = 168 ; + typeOfSecondFixedSurface = 255 ; + scaledValueOfSecondFixedSurface = missing() ; + scaleFactorOfSecondFixedSurface = missing() ; typeOfStatisticalProcessing = 0 ; } #Time-mean sea water potential temperature tendency due to newtonian relaxation diff --git a/definitions/grib2/cfVarName.legacy.def b/definitions/grib2/cfVarName.legacy.def index 40534b276..38c2868e6 100644 --- a/definitions/grib2/cfVarName.legacy.def +++ b/definitions/grib2/cfVarName.legacy.def @@ -157,3 +157,20 @@ typeOfStatisticalProcessing = 1 ; lengthOfTimeRange = 24 ; } +#Upward sea water velocity +'wo' = { + discipline = 10 ; + parameterCategory = 4 ; + parameterNumber = 27 ; + typeOfFirstFixedSurface = 168 ; + typeOfSecondFixedSurface = 168 ; +} +#Time-mean upward sea water velocity +'avg_wo' = { + discipline = 10 ; + parameterCategory = 4 ; + parameterNumber = 27 ; + typeOfFirstFixedSurface = 168 ; + typeOfSecondFixedSurface = 168 ; + typeOfStatisticalProcessing = 0 ; +} diff --git a/definitions/grib2/localConcepts/s2s/cfVarName.def b/definitions/grib2/localConcepts/s2s/cfVarName.def index 13a346191..156261d3c 100644 --- a/definitions/grib2/localConcepts/s2s/cfVarName.def +++ b/definitions/grib2/localConcepts/s2s/cfVarName.def @@ -182,4 +182,11 @@ scaledValueOfSecondFixedSurface = missing() ; scaleFactorOfSecondFixedSurface = missing() ; typeOfStatisticalProcessing = 0 ; + } +#Snow depth water equivalent +'sd' = { + discipline = 0 ; + parameterCategory = 1 ; + parameterNumber = 60 ; + typeOfStatisticalProcessing = 0 ; } diff --git a/definitions/grib2/localConcepts/s2s/name.def b/definitions/grib2/localConcepts/s2s/name.def index 67dd7608e..cc091e047 100644 --- a/definitions/grib2/localConcepts/s2s/name.def +++ b/definitions/grib2/localConcepts/s2s/name.def @@ -182,4 +182,11 @@ scaledValueOfSecondFixedSurface = missing() ; scaleFactorOfSecondFixedSurface = missing() ; typeOfStatisticalProcessing = 0 ; + } +#Snow depth water equivalent +'Snow depth water equivalent' = { + discipline = 0 ; + parameterCategory = 1 ; + parameterNumber = 60 ; + typeOfStatisticalProcessing = 0 ; } diff --git a/definitions/grib2/localConcepts/s2s/paramId.def b/definitions/grib2/localConcepts/s2s/paramId.def index 1156927c0..8db27e06d 100644 --- a/definitions/grib2/localConcepts/s2s/paramId.def +++ b/definitions/grib2/localConcepts/s2s/paramId.def @@ -182,4 +182,11 @@ scaledValueOfSecondFixedSurface = missing() ; scaleFactorOfSecondFixedSurface = missing() ; typeOfStatisticalProcessing = 0 ; + } +#Snow depth water equivalent +'228141' = { + discipline = 0 ; + parameterCategory = 1 ; + parameterNumber = 60 ; + typeOfStatisticalProcessing = 0 ; } diff --git a/definitions/grib2/localConcepts/s2s/shortName.def b/definitions/grib2/localConcepts/s2s/shortName.def index 824ef22e7..025f7fe6f 100644 --- a/definitions/grib2/localConcepts/s2s/shortName.def +++ b/definitions/grib2/localConcepts/s2s/shortName.def @@ -182,4 +182,11 @@ scaledValueOfSecondFixedSurface = missing() ; scaleFactorOfSecondFixedSurface = missing() ; typeOfStatisticalProcessing = 0 ; + } +#Snow depth water equivalent +'sd' = { + discipline = 0 ; + parameterCategory = 1 ; + parameterNumber = 60 ; + typeOfStatisticalProcessing = 0 ; } diff --git a/definitions/grib2/localConcepts/s2s/units.def b/definitions/grib2/localConcepts/s2s/units.def index 4b7e30d55..17ff29efa 100644 --- a/definitions/grib2/localConcepts/s2s/units.def +++ b/definitions/grib2/localConcepts/s2s/units.def @@ -182,4 +182,11 @@ scaledValueOfSecondFixedSurface = missing() ; scaleFactorOfSecondFixedSurface = missing() ; typeOfStatisticalProcessing = 0 ; + } +#Snow depth water equivalent +'kg m**-2' = { + discipline = 0 ; + parameterCategory = 1 ; + parameterNumber = 60 ; + typeOfStatisticalProcessing = 0 ; } diff --git a/definitions/grib2/name.def b/definitions/grib2/name.def index 2a8c0a64b..b3b91952a 100644 --- a/definitions/grib2/name.def +++ b/definitions/grib2/name.def @@ -7841,7 +7841,9 @@ parameterCategory = 4 ; parameterNumber = 27 ; typeOfFirstFixedSurface = 168 ; - typeOfSecondFixedSurface = 168 ; + typeOfSecondFixedSurface = 255 ; + scaledValueOfSecondFixedSurface = missing() ; + scaleFactorOfSecondFixedSurface = missing() ; } #Sea water potential temperature tendency due to newtonian relaxation 'Sea water potential temperature tendency due to newtonian relaxation' = { @@ -9028,7 +9030,9 @@ parameterCategory = 4 ; parameterNumber = 27 ; typeOfFirstFixedSurface = 168 ; - typeOfSecondFixedSurface = 168 ; + typeOfSecondFixedSurface = 255 ; + scaledValueOfSecondFixedSurface = missing() ; + scaleFactorOfSecondFixedSurface = missing() ; typeOfStatisticalProcessing = 0 ; } #Time-mean sea water potential temperature tendency due to newtonian relaxation diff --git a/definitions/grib2/name.legacy.def b/definitions/grib2/name.legacy.def index 62e965b1c..d55237eb9 100644 --- a/definitions/grib2/name.legacy.def +++ b/definitions/grib2/name.legacy.def @@ -157,3 +157,20 @@ typeOfStatisticalProcessing = 1 ; lengthOfTimeRange = 24 ; } +#Upward sea water velocity +'Upward sea water velocity' = { + discipline = 10 ; + parameterCategory = 4 ; + parameterNumber = 27 ; + typeOfFirstFixedSurface = 168 ; + typeOfSecondFixedSurface = 168 ; +} +#Time-mean upward sea water velocity +'Time-mean upward sea water velocity' = { + discipline = 10 ; + parameterCategory = 4 ; + parameterNumber = 27 ; + typeOfFirstFixedSurface = 168 ; + typeOfSecondFixedSurface = 168 ; + typeOfStatisticalProcessing = 0 ; +} diff --git a/definitions/grib2/paramId.def b/definitions/grib2/paramId.def index 453a51920..4c79c56a5 100644 --- a/definitions/grib2/paramId.def +++ b/definitions/grib2/paramId.def @@ -7841,7 +7841,9 @@ parameterCategory = 4 ; parameterNumber = 27 ; typeOfFirstFixedSurface = 168 ; - typeOfSecondFixedSurface = 168 ; + typeOfSecondFixedSurface = 255 ; + scaledValueOfSecondFixedSurface = missing() ; + scaleFactorOfSecondFixedSurface = missing() ; } #Sea water potential temperature tendency due to newtonian relaxation '262508' = { @@ -9028,7 +9030,9 @@ parameterCategory = 4 ; parameterNumber = 27 ; typeOfFirstFixedSurface = 168 ; - typeOfSecondFixedSurface = 168 ; + typeOfSecondFixedSurface = 255 ; + scaledValueOfSecondFixedSurface = missing() ; + scaleFactorOfSecondFixedSurface = missing() ; typeOfStatisticalProcessing = 0 ; } #Time-mean sea water potential temperature tendency due to newtonian relaxation diff --git a/definitions/grib2/paramId.legacy.def b/definitions/grib2/paramId.legacy.def index 489f33c41..ad0ba7521 100644 --- a/definitions/grib2/paramId.legacy.def +++ b/definitions/grib2/paramId.legacy.def @@ -157,3 +157,20 @@ typeOfStatisticalProcessing = 1 ; lengthOfTimeRange = 24 ; } +#Upward sea water velocity +'262507' = { + discipline = 10 ; + parameterCategory = 4 ; + parameterNumber = 27 ; + typeOfFirstFixedSurface = 168 ; + typeOfSecondFixedSurface = 168 ; +} +#Time-mean upward sea water velocity +'263507' = { + discipline = 10 ; + parameterCategory = 4 ; + parameterNumber = 27 ; + typeOfFirstFixedSurface = 168 ; + typeOfSecondFixedSurface = 168 ; + typeOfStatisticalProcessing = 0 ; +} diff --git a/definitions/grib2/shortName.def b/definitions/grib2/shortName.def index 806350301..614a8454d 100644 --- a/definitions/grib2/shortName.def +++ b/definitions/grib2/shortName.def @@ -7841,7 +7841,9 @@ parameterCategory = 4 ; parameterNumber = 27 ; typeOfFirstFixedSurface = 168 ; - typeOfSecondFixedSurface = 168 ; + typeOfSecondFixedSurface = 255 ; + scaledValueOfSecondFixedSurface = missing() ; + scaleFactorOfSecondFixedSurface = missing() ; } #Sea water potential temperature tendency due to newtonian relaxation 'thetaodmp' = { @@ -9028,7 +9030,9 @@ parameterCategory = 4 ; parameterNumber = 27 ; typeOfFirstFixedSurface = 168 ; - typeOfSecondFixedSurface = 168 ; + typeOfSecondFixedSurface = 255 ; + scaledValueOfSecondFixedSurface = missing() ; + scaleFactorOfSecondFixedSurface = missing() ; typeOfStatisticalProcessing = 0 ; } #Time-mean sea water potential temperature tendency due to newtonian relaxation diff --git a/definitions/grib2/shortName.legacy.def b/definitions/grib2/shortName.legacy.def index 40534b276..38c2868e6 100644 --- a/definitions/grib2/shortName.legacy.def +++ b/definitions/grib2/shortName.legacy.def @@ -157,3 +157,20 @@ typeOfStatisticalProcessing = 1 ; lengthOfTimeRange = 24 ; } +#Upward sea water velocity +'wo' = { + discipline = 10 ; + parameterCategory = 4 ; + parameterNumber = 27 ; + typeOfFirstFixedSurface = 168 ; + typeOfSecondFixedSurface = 168 ; +} +#Time-mean upward sea water velocity +'avg_wo' = { + discipline = 10 ; + parameterCategory = 4 ; + parameterNumber = 27 ; + typeOfFirstFixedSurface = 168 ; + typeOfSecondFixedSurface = 168 ; + typeOfStatisticalProcessing = 0 ; +} diff --git a/definitions/grib2/stepUnits.def b/definitions/grib2/stepUnits.def index cffc236be..d1beeb517 100644 --- a/definitions/grib2/stepUnits.def +++ b/definitions/grib2/stepUnits.def @@ -4,7 +4,8 @@ # template_nofail default_step_units "grib2/localConcepts/[centre:s]/default_step_units.def"; # codetable[1] stepUnits 'stepUnits.table' = defaultStepUnits : transient,dump,no_copy; -meta stepUnits optimal_step_units(forecastTime,indicatorOfUnitOfTimeRange,lengthOfTimeRange,indicatorOfUnitForTimeRange) : transient,dump; +# See ECC-1768 re why no_copy is needed +meta stepUnits optimal_step_units(forecastTime,indicatorOfUnitOfTimeRange,lengthOfTimeRange,indicatorOfUnitForTimeRange) : dump,no_copy; transient startStepUnit = 255 : hidden; # 255 means MISSING. See code table 4.4 transient endStepUnit = 255 : hidden; # The lowercase version is to unify it with the helper key in the MARS language diff --git a/definitions/grib2/tables/local/ecmf/1/4.2.0.1.table b/definitions/grib2/tables/local/ecmf/1/4.2.0.1.table index 5e02901a1..1ee04a741 100644 --- a/definitions/grib2/tables/local/ecmf/1/4.2.0.1.table +++ b/definitions/grib2/tables/local/ecmf/1/4.2.0.1.table @@ -1,7 +1,7 @@ # Code table 4.2 - discipline=0 category=1 for ECMWF 192 192 Snow evaporation rate (kg m-2 s-1) 193 193 Total precipitation rate (m s-1) -194 194 Accumulated freezing rain (m) +194 194 Freezing rain precipitation rate (m s-1) 195 195 Convective precipitation rate (m s-1) 196 196 Large-scale precipitation rate (m s-1) 197 197 Snow evaporation rate (m of water equivalent s-1) diff --git a/definitions/grib2/units.def b/definitions/grib2/units.def index 411ec5996..d96d6bd9e 100644 --- a/definitions/grib2/units.def +++ b/definitions/grib2/units.def @@ -7841,7 +7841,9 @@ parameterCategory = 4 ; parameterNumber = 27 ; typeOfFirstFixedSurface = 168 ; - typeOfSecondFixedSurface = 168 ; + typeOfSecondFixedSurface = 255 ; + scaledValueOfSecondFixedSurface = missing() ; + scaleFactorOfSecondFixedSurface = missing() ; } #Sea water potential temperature tendency due to newtonian relaxation 'K s**-1' = { @@ -9028,7 +9030,9 @@ parameterCategory = 4 ; parameterNumber = 27 ; typeOfFirstFixedSurface = 168 ; - typeOfSecondFixedSurface = 168 ; + typeOfSecondFixedSurface = 255 ; + scaledValueOfSecondFixedSurface = missing() ; + scaleFactorOfSecondFixedSurface = missing() ; typeOfStatisticalProcessing = 0 ; } #Time-mean sea water potential temperature tendency due to newtonian relaxation diff --git a/definitions/grib2/units.legacy.def b/definitions/grib2/units.legacy.def index 145f4aaac..f05dae694 100644 --- a/definitions/grib2/units.legacy.def +++ b/definitions/grib2/units.legacy.def @@ -157,3 +157,20 @@ typeOfStatisticalProcessing = 1 ; lengthOfTimeRange = 24 ; } +#Upward sea water velocity +'m s**-1' = { + discipline = 10 ; + parameterCategory = 4 ; + parameterNumber = 27 ; + typeOfFirstFixedSurface = 168 ; + typeOfSecondFixedSurface = 168 ; +} +#Time-mean upward sea water velocity +'m s**-1' = { + discipline = 10 ; + parameterCategory = 4 ; + parameterNumber = 27 ; + typeOfFirstFixedSurface = 168 ; + typeOfSecondFixedSurface = 168 ; + typeOfStatisticalProcessing = 0 ; +} diff --git a/definitions/mars/grib.mmsf.fc.def b/definitions/mars/grib.mmsf.fc.def index d213cee45..8791cb086 100644 --- a/definitions/mars/grib.mmsf.fc.def +++ b/definitions/mars/grib.mmsf.fc.def @@ -1,8 +1,14 @@ -alias mars.step = endStep; +if (levtype is "o2d" || levtype is "o3d") { + alias mars.step = stepRange; +} else { + alias mars.step = endStep; +} + if (class is "od") { alias mars.system = systemNumber; } if (class is "me") { alias mars.system = systemNumber; } if (class is "en") { alias mars.system = systemNumber; } if (class is "c3") { alias mars.system = systemNumber; } +if (class is "ci") { alias mars.system = systemNumber; } alias mars.number = perturbationNumber; alias mars.method = methodNumber; diff --git a/definitions/mars/grib.msmm.fcmean.def b/definitions/mars/grib.msmm.fcmean.def index 800a49134..a13ca9918 100644 --- a/definitions/mars/grib.msmm.fcmean.def +++ b/definitions/mars/grib.msmm.fcmean.def @@ -10,6 +10,7 @@ if (class is "od") { alias mars.system = systemNumber; } if (class is "me") { alias mars.system = systemNumber; } if (class is "en") { alias mars.system = systemNumber; } if (class is "c3") { alias mars.system = systemNumber; } +if (class is "ci") { alias mars.system = systemNumber; } # See ECC-624 if (centre == 80 && subCentre == 98 && class is "c3") { diff --git a/src/grib_iterator_class_lambert_conformal.cc b/src/grib_iterator_class_lambert_conformal.cc index 3e9d51055..5d15a67b6 100644 --- a/src/grib_iterator_class_lambert_conformal.cc +++ b/src/grib_iterator_class_lambert_conformal.cc @@ -85,21 +85,21 @@ static void init_class(grib_iterator_class* c) #define EPSILON 1.0e-10 #ifndef M_PI -#define M_PI 3.14159265358979323846 /* Whole pie */ +#define M_PI 3.14159265358979323846 // Whole pie #endif #ifndef M_PI_2 -#define M_PI_2 1.57079632679489661923 /* Half a pie */ +#define M_PI_2 1.57079632679489661923 // Half a pie #endif #ifndef M_PI_4 -#define M_PI_4 0.78539816339744830962 /* Quarter of a pie */ +#define M_PI_4 0.78539816339744830962 // Quarter of a pie #endif -#define RAD2DEG 57.29577951308232087684 /* 180 over pi */ -#define DEG2RAD 0.01745329251994329576 /* pi over 180 */ +#define RAD2DEG 57.29577951308232087684 // 180 over pi +#define DEG2RAD 0.01745329251994329576 // pi over 180 -/* Adjust longitude (in radians) to range -180 to 180 */ +// Adjust longitude (in radians) to range -180 to 180 static double adjust_lon_radians(double lon) { if (lon > M_PI) lon -= 2 * M_PI; @@ -107,16 +107,15 @@ static double adjust_lon_radians(double lon) return lon; } -/* Function to compute the latitude angle, phi2, for the inverse - * From the book "Map Projections-A Working Manual-John P. Snyder (1987)" - * Equation (7-9) involves rapidly converging iteration: Calculate t from (15-11) - * Then, assuming an initial trial phi equal to (pi/2 - 2*arctan t) in the right side of equation (7-9), - * calculate phi on the left side. Substitute the calculated phi) into the right side, - * calculate a new phi, etc., until phi does not change significantly from the preceding trial value of phi - */ +// Function to compute the latitude angle, phi2, for the inverse +// From the book "Map Projections-A Working Manual-John P. Snyder (1987)" +// Equation (7-9) involves rapidly converging iteration: Calculate t from (15-11) +// Then, assuming an initial trial phi equal to (pi/2 - 2*arctan t) in the right side of equation (7-9), +// calculate phi on the left side. Substitute the calculated phi) into the right side, +// calculate a new phi, etc., until phi does not change significantly from the preceding trial value of phi static double compute_phi( - double eccent, /* Spheroid eccentricity */ - double ts, /* Constant value t */ + double eccent, // Spheroid eccentricity + double ts, // Constant value t int* error) { double eccnth, phi, con, dphi, sinpi; @@ -136,19 +135,19 @@ static double compute_phi( return 0; } -/* Compute the constant small m which is the radius of - a parallel of latitude, phi, divided by the semimajor axis */ +// Compute the constant small m which is the radius of +// a parallel of latitude, phi, divided by the semimajor axis static double compute_m(double eccent, double sinphi, double cosphi) { const double con = eccent * sinphi; return ((cosphi / (sqrt(1.0 - con * con)))); } -/* Compute the constant small t for use in the forward computations */ +// Compute the constant small t for use in the forward computations static double compute_t( - double eccent, /* Eccentricity of the spheroid */ - double phi, /* Latitude phi */ - double sinphi) /* Sine of the latitude */ + double eccent, // Eccentricity of the spheroid + double phi, // Latitude phi + double sinphi) // Sine of the latitude { double con = eccent * sinphi; double com = 0.5 * eccent; @@ -162,7 +161,34 @@ static double calculate_eccentricity(double minor, double major) return sqrt(1.0 - temp * temp); } -static int init_sphere(grib_handle* h, +static void xy2lonlat(double radius, double n, double f, double rho0_bare, double LoVInRadians, + double x, double y, + double* lonDeg, double* latDeg) +{ + DEBUG_ASSERT(radius > 0); + DEBUG_ASSERT(n != 0.0); + x /= radius; + y /= radius; + y = rho0_bare - y; + double rho = hypot(x, y); + if (rho != 0.0) { + if (n < 0.0) { + rho = -rho; + x = -x; + y = -y; + } + double latRadians = 2. * atan(pow(f / rho, 1.0/n)) - M_PI_2; + double lonRadians = atan2(x, y) / n; + *lonDeg = (lonRadians + LoVInRadians) * RAD2DEG; + *latDeg = latRadians * RAD2DEG; + } + else { + *lonDeg = 0.0; + *latDeg = (n > 0.0 ? M_PI_2 : -M_PI_2) * RAD2DEG; + } +} + +static int init_sphere(const grib_handle* h, grib_iterator_lambert_conformal* self, size_t nv, long nx, long ny, double LoVInDegrees, @@ -171,9 +197,7 @@ static int init_sphere(grib_handle* h, double LoVInRadians, double Latin1InRadians, double Latin2InRadians, double LaDInRadians) { - int i, j; - double f, n, rho, rho0, angle, x0, y0, x, y, tmp, tmp2; - double latDeg, lonDeg, lonDiff; + double n, x, y; if (fabs(Latin1InRadians - Latin2InRadians) < 1E-09) { n = sin(Latin1InRadians); @@ -182,26 +206,25 @@ static int init_sphere(grib_handle* h, log(tan(M_PI_4 + Latin2InRadians / 2.0) / tan(M_PI_4 + Latin1InRadians / 2.0)); } - f = (cos(Latin1InRadians) * pow(tan(M_PI_4 + Latin1InRadians / 2.0), n)) / n; - rho = radius * f * pow(tan(M_PI_4 + latFirstInRadians / 2.0), -n); - rho0 = radius * f * pow(tan(M_PI_4 + LaDInRadians / 2.0), -n); - if (n < 0) /* adjustment for southern hemisphere */ - rho0 = -rho0; - lonDiff = lonFirstInRadians - LoVInRadians; + double f = (cos(Latin1InRadians) * pow(tan(M_PI_4 + Latin1InRadians / 2.0), n)) / n; + double rho = radius * f * pow(tan(M_PI_4 + latFirstInRadians / 2.0), -n); + double rho0_bare = f * pow(tan(M_PI_4 + LaDInRadians / 2.0), -n); + double rho0 = radius * rho0_bare; // scaled + double lonDiff = lonFirstInRadians - LoVInRadians; - /* Adjust longitude to range -180 to 180 */ + // Adjust longitude to range -180 to 180 if (lonDiff > M_PI) lonDiff -= 2 * M_PI; if (lonDiff < -M_PI) lonDiff += 2 * M_PI; - angle = n * lonDiff; - x0 = rho * sin(angle); - y0 = rho0 - rho * cos(angle); - /*Dx = iScansNegatively == 0 ? Dx : -Dx;*/ - /* GRIB-405: Don't change sign of Dy. Latitudes ALWAYS increase from latitudeOfFirstGridPoint */ - /*Dy = jScansPositively == 1 ? Dy : -Dy;*/ - - /* Allocate latitude and longitude arrays */ + double angle = n * lonDiff; + double x0 = rho * sin(angle); + double y0 = rho0 - rho * cos(angle); + // Dx = iScansNegatively == 0 ? Dx : -Dx; + // GRIB-405: Don't change sign of Dy. Latitudes ALWAYS increase from latitudeOfFirstGridPoint + // Dy = jScansPositively == 1 ? Dy : -Dy; + + // Allocate latitude and longitude arrays self->lats = (double*)grib_context_malloc(h->context, nv * sizeof(double)); if (!self->lats) { grib_context_log(h->context, GRIB_LOG_ERROR, "%s: Error allocating %zu bytes", ITER, nv * sizeof(double)); @@ -213,20 +236,32 @@ static int init_sphere(grib_handle* h, return GRIB_OUT_OF_MEMORY; } - /* Populate our arrays */ - for (j = 0; j < ny; j++) { + double latDeg = 0, lonDeg = 0; + for (long j = 0; j < ny; j++) { y = y0 + j * Dy; - if (n < 0) { /* adjustment for southern hemisphere */ - y = -y; + for (long i = 0; i < nx; i++) { + const long index = i + j * nx; + x = x0 + i * Dx; + xy2lonlat(radius, n, f, rho0_bare, LoVInRadians, x, y, &lonDeg, &latDeg); + self->lons[index] = lonDeg; + self->lats[index] = latDeg; } + } + +#if 0 + for (j = 0; j < ny; j++) { + y = y0 + j * Dy; + //if (n < 0) { /* adjustment for southern hemisphere */ + // y = -y; + //} tmp = rho0 - y; tmp2 = tmp * tmp; for (i = 0; i < nx; i++) { int index = i + j * nx; - x = x0 + i * Dx; - if (n < 0) { /* adjustment for southern hemisphere */ - x = -x; - } + x = x0 + i * Dx; + //if (n < 0) { /* adjustment for southern hemisphere */ + // x = -x; + //} angle = atan2(x, tmp); /* See ECC-524 */ rho = sqrt(x * x + tmp2); if (n <= 0) rho = -rho; @@ -237,12 +272,12 @@ static int init_sphere(grib_handle* h, self->lats[index] = latDeg; } } - +#endif return GRIB_SUCCESS; } -/* Oblate spheroid */ -static int init_oblate(grib_handle* h, +// Oblate spheroid +static int init_oblate(const grib_handle* h, grib_iterator_lambert_conformal* self, size_t nv, long nx, long ny, double LoVInDegrees, @@ -254,20 +289,20 @@ static int init_oblate(grib_handle* h, { int i, j, err = 0; double x0, y0, x, y, latRad, lonRad, latDeg, lonDeg, sinphi, ts, rh1, theta; - double false_easting; /* x offset in meters */ - double false_northing; /* y offset in meters */ - - double ns; /* ratio of angle between meridian */ - double F; /* flattening of ellipsoid */ - double rh; /* height above ellipsoid */ - double sin_po; /* sin value */ - double cos_po; /* cos value */ - double con; /* temporary variable */ - double ms1; /* small m 1 */ - double ms2; /* small m 2 */ - double ts0; /* small t 0 */ - double ts1; /* small t 1 */ - double ts2; /* small t 2 */ + double false_easting; // x offset in meters + double false_northing; // y offset in meters + + double ns; // ratio of angle between meridian + double F; // flattening of ellipsoid + double rh; // height above ellipsoid + double sin_po; // sin value + double cos_po; // cos value + double con; // temporary variable + double ms1; // small m 1 + double ms2; // small m 2 + double ts0; // small t 0 + double ts1; // small t 1 + double ts2; // small t 2 double e = calculate_eccentricity(earthMinorAxisInMetres, earthMajorAxisInMetres); @@ -292,7 +327,7 @@ static int init_oblate(grib_handle* h, F = ms1 / (ns * pow(ts1, ns)); rh = earthMajorAxisInMetres * F * pow(ts0, ns); - /* Forward projection: convert lat,lon to x,y */ + // Forward projection: convert lat,lon to x,y con = fabs(fabs(latFirstInRadians) - M_PI_2); if (con > EPSILON) { sinphi = sin(latFirstInRadians); @@ -312,7 +347,7 @@ static int init_oblate(grib_handle* h, x0 = -x0; y0 = -y0; - /* Allocate latitude and longitude arrays */ + // Allocate latitude and longitude arrays self->lats = (double*)grib_context_malloc(h->context, nv * sizeof(double)); if (!self->lats) { grib_context_log(h->context, GRIB_LOG_ERROR, "%s: Error allocating %zu bytes", ITER, nv * sizeof(double)); @@ -324,7 +359,7 @@ static int init_oblate(grib_handle* h, return GRIB_OUT_OF_MEMORY; } - /* Populate our arrays */ + // Populate our arrays false_easting = x0; false_northing = y0; for (j = 0; j < ny; j++) { @@ -333,7 +368,7 @@ static int init_oblate(grib_handle* h, const int index = i + j * nx; double _x, _y; x = i * Dx; - /* Inverse projection to convert from x,y to lat,lon */ + // Inverse projection to convert from x,y to lat,lon _x = x - false_easting; _y = rh - y + false_northing; rh1 = sqrt(_x * _x + _y * _y); @@ -363,7 +398,7 @@ static int init_oblate(grib_handle* h, if (i == 0 && j == 0) { DEBUG_ASSERT(fabs(latFirstInRadians - latRad) <= EPSILON); } - latDeg = latRad * RAD2DEG; /* Convert to degrees */ + latDeg = latRad * RAD2DEG; // Convert to degrees lonDeg = normalise_longitude_in_degrees(lonRad * RAD2DEG); self->lons[index] = lonDeg; self->lats[index] = latDeg; @@ -393,7 +428,7 @@ static int init(grib_iterator* iter, grib_handle* h, grib_arguments* args) const char* sLatin2InDegrees = grib_arguments_get_name(h, args, self->carg++); const char* slatFirstInDegrees = grib_arguments_get_name(h, args, self->carg++); const char* slonFirstInDegrees = grib_arguments_get_name(h, args, self->carg++); - /* Dx and Dy are in Metres */ + // Dx and Dy are in Metres const char* sDx = grib_arguments_get_name(h, args, self->carg++); const char* sDy = grib_arguments_get_name(h, args, self->carg++); const char* siScansNegatively = grib_arguments_get_name(h, args, self->carg++); @@ -443,16 +478,16 @@ static int init(grib_iterator* iter, grib_handle* h, grib_arguments* args) if ((err = grib_get_long_internal(h, salternativeRowScanning, &alternativeRowScanning)) != GRIB_SUCCESS) return err; - /* Standard Parallels cannot be equal and on opposite sides of the equator */ + // Standard Parallels cannot be equal and on opposite sides of the equator if (fabs(Latin1InDegrees + Latin2InDegrees) < EPSILON) { grib_context_log(h->context, GRIB_LOG_ERROR, "%s: Cannot have equal latitudes for standard parallels on opposite sides of equator", ITER); return GRIB_WRONG_GRID; } - /* - * See Wolfram MathWorld: http://mathworld.wolfram.com/LambertConformalConicProjection.html - */ + // + // See Wolfram MathWorld: http://mathworld.wolfram.com/LambertConformalConicProjection.html + // latFirstInRadians = latFirstInDegrees * DEG2RAD; lonFirstInRadians = lonFirstInDegrees * DEG2RAD; Latin1InRadians = Latin1InDegrees * DEG2RAD; @@ -478,7 +513,7 @@ static int init(grib_iterator* iter, grib_handle* h, grib_arguments* args) iter->e = -1; - /* Apply the scanning mode flags which may require data array to be transformed */ + // Apply the scanning mode flags which may require data array to be transformed err = transform_iterator_data(h->context, iter->data, iScansNegatively, jScansPositively, jPointsAreConsecutive, alternativeRowScanning, iter->nv, nx, ny); diff --git a/src/step_unit.cc b/src/step_unit.cc index bd192a6e7..a842430f2 100644 --- a/src/step_unit.cc +++ b/src/step_unit.cc @@ -12,8 +12,6 @@ namespace eccodes { -Unit::Map Unit::map_{}; - std::vector Unit::grib_selected_units = { Unit::Value::SECOND, Unit::Value::MINUTE, @@ -39,7 +37,7 @@ std::vector Unit::complete_unit_order_ = { }; template <> long Unit::value() const { - return map_.unit_to_long(internal_value_); + return get_converter().unit_to_long(internal_value_); } template <> Unit::Value Unit::value() const { @@ -47,7 +45,7 @@ template <> Unit::Value Unit::value() const { } template <> std::string Unit::value() const { - return map_.unit_to_name(internal_value_); + return get_converter().unit_to_name(internal_value_); } } // namespace eccodes diff --git a/src/step_unit.h b/src/step_unit.h index bd078f403..0d6bc4e4e 100644 --- a/src/step_unit.h +++ b/src/step_unit.h @@ -68,7 +68,7 @@ class Unit { explicit Unit(const std::string& unit_value) { try { - internal_value_ = map_.name_to_unit(unit_value); + internal_value_ = get_converter().name_to_unit(unit_value); } catch (std::exception& e) { throw std::runtime_error(std::string{"Unit not found "} + e.what()); } @@ -76,15 +76,15 @@ class Unit { explicit Unit(long unit_value) { try { - internal_value_ = map_.long_to_unit(unit_value); + internal_value_ = get_converter().long_to_unit(unit_value); } catch (std::exception& e) { throw std::runtime_error(std::string{"Unit not found "} + e.what()); } } - bool operator>(const Unit& other) const {return map_.unit_to_duration(internal_value_) > map_.unit_to_duration(other.internal_value_);} - bool operator==(const Value value) const {return map_.unit_to_duration(internal_value_) == map_.unit_to_duration(value);} - bool operator==(const Unit& unit) const {return map_.unit_to_duration(internal_value_) == map_.unit_to_duration(unit.internal_value_);} + bool operator>(const Unit& other) const {return get_converter().unit_to_duration(internal_value_) > get_converter().unit_to_duration(other.internal_value_);} + bool operator==(const Value value) const {return get_converter().unit_to_duration(internal_value_) == get_converter().unit_to_duration(value);} + bool operator==(const Unit& unit) const {return get_converter().unit_to_duration(internal_value_) == get_converter().unit_to_duration(unit.internal_value_);} bool operator!=(const Unit& unit) const {return !(*this == unit);} bool operator!=(const Value value) const {return !(*this == value);} @@ -177,13 +177,14 @@ class Unit { Value internal_value_; - static Map map_; public: - static Map& get_converter() {return map_;} + static Map& get_converter() { + static Map map_; + return map_; + } }; - template Seconds to_seconds(long value, const Unit& unit) { Seconds seconds; diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt index 71babe949..32c20a5a5 100644 --- a/tests/CMakeLists.txt +++ b/tests/CMakeLists.txt @@ -254,6 +254,7 @@ if( HAVE_BUILD_TOOLS ) grib_ecc-1000 grib_ecc-1001 grib_ecc-1030 + grib_ecc-1364 grib_ecc-1397 grib_ecc-1425 grib_ecc-1467 diff --git a/tests/bufr_compare.sh b/tests/bufr_compare.sh index 882fa3293..27f0c071b 100755 --- a/tests/bufr_compare.sh +++ b/tests/bufr_compare.sh @@ -351,6 +351,31 @@ set -e [ $status -ne 0 ] grep -q "unreadable message" $fLog +###??? +f1="aaen_55.bufr" +f2="aaen_55.bufr" +set +e +${tools_dir}/bufr_compare -H -c edition $f1 $f2 > $fLog 2>&1 +status=$? +set -e +[ $status -ne 0 ] +grep -q "options are incompatible" $fLog + +set +e +${tools_dir}/bufr_compare -a edition $f1 $f2 > $fLog 2>&1 +status=$? +set -e +[ $status -ne 0 ] +grep -q "a option requires -c option" $fLog + + +set +e +${tools_dir}/bufr_compare nosuchfile $f1 > $fLog 2>&1 +status=$? +set -e +[ $status -ne 0 ] + + # Clean up # ------------- diff --git a/tests/grib_ecc-1364.sh b/tests/grib_ecc-1364.sh new file mode 100755 index 000000000..cb4555834 --- /dev/null +++ b/tests/grib_ecc-1364.sh @@ -0,0 +1,52 @@ +#!/bin/sh +# (C) Copyright 2005- 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. +# +# 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. +# + +. ./include.ctest.sh + + +# ECC-1364: GRIB: Geoiterator for Lambert Conformal in the southern hemisphere + +label="grib_ecc-1364_test" +tempGrib=temp.$label.grib +tempFilt=temp.$label.filt +tempLog=temp.$label.log +tempRef=temp.$label.ref + +sample_grib1=$ECCODES_SAMPLES_PATH/GRIB1.tmpl + +# Create a GRIB with a similar grid to the one in the JIRA issue +cat >$tempFilt< $tempLog + +${tools_dir}/grib_ls -l -11.6277,-47.9583,1 $tempGrib > $tempLog +grep -q "Grid Point chosen #1 index=1247750 " $tempLog +grep -q "index=1247750 .* distance=0.01 " $tempLog + + +# Clean up +rm -f $tempGrib $tempFilt $tempLog $tempRef diff --git a/tests/grib_sub_hourly.sh b/tests/grib_sub_hourly.sh index 75b113b7e..fdd174c28 100755 --- a/tests/grib_sub_hourly.sh +++ b/tests/grib_sub_hourly.sh @@ -520,6 +520,20 @@ cat $tempFilt ${tools_dir}/grib_filter $tempFilt $data_dir/constant_field.grib2 unset ECCODES_GRIB_HOURLY_STEPS_WITH_UNITS + +# Changing the product definition template +# ---------------------------------------- +# See ECC-1768 +${tools_dir}/grib_set -s step=62m $sample_g2 $temp +${tools_dir}/grib_set -s productDefinitionTemplateNumber=8 $temp $temp2 + +${tools_dir}/grib_set -s productDefinitionTemplateNumber=8,stepUnits=s,step=0 $sample_g2 $temp +grib_check_key_equals $temp '-p stepUnits:s,startStep,productDefinitionTemplateNumber' 's 0s 8' + +${tools_dir}/grib_set -s productDefinitionTemplateNumber=8,stepUnits=m,step=60 $sample_g2 $temp +grib_check_key_equals $temp '-p stepUnits:s,productDefinitionTemplateNumber' 'h 8' + + # Bad stepUnits set +e ${tools_dir}/grib_set -s stepUnits=190 $sample_g2 $temp > $tempText 2>&1