Commit 4f86c67c authored by Uwe Schulzweida's avatar Uwe Schulzweida
Browse files

gridInqUvRelativeToGrid(): added cgribex support.

parent d43badb1
......@@ -2,6 +2,10 @@
* Version 1.8.1 released
2017-02-27 Uwe Schulzweida
* gridInqUvRelativeToGrid(): added cgribex support
2017-02-24 Uwe Schulzweida
* Added function gridInqUvRelativeToGrid() (patch from Michal Koutek, KMNI)
......
......@@ -713,10 +713,8 @@ void gribapiGetGrid(grib_handle *gh, grid_t *grid)
Error("Unsupported grid type: %s", gridNamePtr(gridtype));
}
if ( gridtype == GRID_LONLAT || gridtype == GRID_CURVILINEAR || gridtype == GRID_LCC )
{
GRIB_CHECK(grib_get_long(gh, "uvRelativeToGrid", &grid->uvRelativeToGrid), 0);
}
if ( gridtype == GRID_GAUSSIAN || gridtype == GRID_LONLAT || gridtype == GRID_PROJECTION || gridtype == GRID_LCC )
GRIB_CHECK(grib_get_long(gh, "uvRelativeToGrid", &grid->uvRelativeToGrid), 0);
grid->type = gridtype;
grid->projtype = projtype;
......
......@@ -29,6 +29,8 @@ typedef struct {
#if defined (HAVE_LIBCGRIBEX)
static bool bit_get(int number, int bit) { return (bool)((number >> bit) & 1); }
static
int cgribexGetGridType(int *isec2)
{
......@@ -163,6 +165,9 @@ void cgribexGetGrid(stream_t *streamptr, int *isec2, int *isec4, grid_t *grid, i
if ( gridtype == GRID_LONLAT || gridtype == GRID_GAUSSIAN || projtype == CDI_PROJ_RLL )
{
bool ijDirectionIncrementGiven = bit_get(ISEC2_ResFlag, 7);
bool uvRelativeToGrid = bit_get(ISEC2_ResFlag, 3);
if ( uvRelativeToGrid ) grid->uvRelativeToGrid = 1;
if ( ISEC4_NumValues != ISEC2_NumLon*ISEC2_NumLat )
Error("numberOfPoints (%d) and gridSize (%d) differ!", ISEC4_NumValues, ISEC2_NumLon*ISEC2_NumLat);
grid->size = ISEC4_NumValues;
......@@ -180,7 +185,7 @@ void cgribexGetGrid(stream_t *streamptr, int *isec2, int *isec4, grid_t *grid, i
if ( ISEC2_LastLon < ISEC2_FirstLon && ISEC2_LastLon < 0 ) ISEC2_LastLon += 360000;
if ( ISEC2_ResFlag && ISEC2_LonIncr > 0 )
if ( ijDirectionIncrementGiven && ISEC2_LonIncr > 0 )
{
if ( abs(ISEC2_LastLon - (ISEC2_FirstLon+ISEC2_LonIncr*(grid->x.size-1))) <= 2 )
{
......@@ -213,7 +218,7 @@ void cgribexGetGrid(stream_t *streamptr, int *isec2, int *isec4, grid_t *grid, i
if ( grid->y.size > 1 && compyinc )
{
bool recompinc = true;
if ( ISEC2_ResFlag && ISEC2_LatIncr > 0 )
if ( ijDirectionIncrementGiven && ISEC2_LatIncr > 0 )
{
if ( abs(ISEC2_LastLat - (ISEC2_FirstLat+ISEC2_LatIncr*(grid->y.size-1))) <= 2 )
{
......@@ -232,6 +237,9 @@ void cgribexGetGrid(stream_t *streamptr, int *isec2, int *isec4, grid_t *grid, i
}
else if ( gridtype == GRID_GAUSSIAN_REDUCED )
{
bool ijDirectionIncrementGiven = bit_get(ISEC2_ResFlag, 7);
bool uvRelativeToGrid = bit_get(ISEC2_ResFlag, 3);
if ( uvRelativeToGrid ) grid->uvRelativeToGrid = 1;
grid->np = ISEC2_NumPar;
grid->size = ISEC4_NumValues;
grid->rowlon = ISEC2_RowLonPtr;
......@@ -246,7 +254,7 @@ void cgribexGetGrid(stream_t *streamptr, int *isec2, int *isec4, grid_t *grid, i
{
if ( ISEC2_LastLon < ISEC2_FirstLon && ISEC2_LastLon < 0 ) ISEC2_LastLon += 360000;
if ( ISEC2_ResFlag && ISEC2_LonIncr > 0 )
if ( ijDirectionIncrementGiven && ISEC2_LonIncr > 0 )
grid->x.inc = ISEC2_LonIncr * 0.001;
else
grid->x.inc = (ISEC2_LastLon - ISEC2_FirstLon) * 0.001 / (grid->x.size - 1);
......@@ -260,7 +268,7 @@ void cgribexGetGrid(stream_t *streamptr, int *isec2, int *isec4, grid_t *grid, i
{
if ( grid->y.size > 1 )
{
if ( ISEC2_ResFlag && ISEC2_LatIncr > 0 )
if ( ijDirectionIncrementGiven && ISEC2_LatIncr > 0 )
grid->y.inc = ISEC2_LatIncr * 0.001;
else
grid->y.inc = (ISEC2_LastLat - ISEC2_FirstLat) * 0.001 / (grid->y.size - 1);
......@@ -272,6 +280,9 @@ void cgribexGetGrid(stream_t *streamptr, int *isec2, int *isec4, grid_t *grid, i
}
else if ( gridtype == GRID_LCC )
{
bool uvRelativeToGrid = bit_get(ISEC2_ResFlag, 3);
if ( uvRelativeToGrid ) grid->uvRelativeToGrid = 1;
if ( ISEC4_NumValues != ISEC2_NumLon*ISEC2_NumLat )
Error("numberOfPoints (%d) and gridSize (%d) differ!",
ISEC4_NumValues, ISEC2_NumLon*ISEC2_NumLat);
......@@ -1719,7 +1730,8 @@ void cgribexDefGrid(int *isec1, int *isec2, double *fsec2, int *isec4, int gridI
if ( ISEC2_NumLon == 1 && ISEC2_NumLat > 1 )
if ( ISEC2_LonIncr == 0 && ISEC2_LatIncr != 0 ) ISEC2_LonIncr = ISEC2_LatIncr;
ISEC2_ResFlag = ( ISEC2_LatIncr == 0 || ISEC2_LonIncr == 0 ) ? 0 : 128;
ISEC2_ResFlag = ( ISEC2_LatIncr == 0 || ISEC2_LonIncr == 0 ) ? 0 : 128; // Set bit 7
if ( gridInqUvRelativeToGrid(gridID) ) ISEC2_ResFlag += 8; // Set bit 3
if ( lrotated )
{
......@@ -1767,6 +1779,9 @@ void cgribexDefGrid(int *isec1, int *isec2, double *fsec2, int *isec4, int gridI
ISEC2_Lambert_ProjFlag = projflag;
ISEC2_ScanFlag = scanflag;
ISEC2_ResFlag = ( ISEC2_Lambert_dx == 0 || ISEC2_Lambert_dy == 0 ) ? 0 : 128; // Set bit 7
if ( gridInqUvRelativeToGrid(gridID) ) ISEC2_ResFlag += 8; // Set bit 3
break;
}
case GRID_SPECTRAL:
......
......@@ -2086,6 +2086,9 @@ void gribapiDefGrid(int editionNumber, grib_handle *gh, int gridID, int comptype
if ( lieee )
GRIB_CHECK(my_grib_set_long(gh, "precision", datatype == CDI_DATATYPE_FLT64 ? 2 : 1), 0);
long uvRelativeToGrid = gridInqUvRelativeToGrid(gridID);
if ( uvRelativeToGrid ) GRIB_CHECK(my_grib_set_long(gh, "uvRelativeToGrid", uvRelativeToGrid), 0);
break;
}
case GRID_LCC:
......@@ -2100,15 +2103,13 @@ void gribapiDefGrid(int editionNumber, grib_handle *gh, int gridID, int comptype
&projflag, &scanflag);
static const char mesg[] = "lambert";
size_t len = sizeof (mesg) -1;
size_t len = sizeof(mesg) -1;
GRIB_CHECK(my_grib_set_string(gh, "gridType", mesg, &len), 0);
GRIB_CHECK(my_grib_set_long(gh, "Nx", xsize), 0);
GRIB_CHECK(my_grib_set_long(gh, "Ny", ysize), 0);
/* FIXME: lround should probably be round here */
GRIB_CHECK(my_grib_set_double(gh, "DxInMetres", (double)lround(xincm)), 0);
/* FIXME: lround should probably be round here */
GRIB_CHECK(my_grib_set_double(gh, "DyInMetres", (double)lround(yincm)), 0);
GRIB_CHECK(my_grib_set_double(gh, "longitudeOfFirstGridPointInDegrees", originLon), 0);
GRIB_CHECK(my_grib_set_double(gh, "latitudeOfFirstGridPointInDegrees", originLat), 0);
......@@ -2122,12 +2123,15 @@ void gribapiDefGrid(int editionNumber, grib_handle *gh, int gridID, int comptype
GRIB_CHECK(my_grib_set_long(gh, "scanningMode", scanflag), 0);
}
long uvRelativeToGrid = gridInqUvRelativeToGrid(gridID);
if ( uvRelativeToGrid ) GRIB_CHECK(my_grib_set_long(gh, "uvRelativeToGrid", uvRelativeToGrid), 0);
break;
}
case GRID_SPECTRAL:
{
static const char mesg[] = "sh";
size_t len = sizeof (mesg) -1;
size_t len = sizeof(mesg) -1;
GRIB_CHECK(my_grib_set_string(gh, "gridType", mesg, &len), 0);
int trunc = gridInqTrunc(gridID);
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment