Commit 3ee56a85 authored by Uwe Schulzweida's avatar Uwe Schulzweida
Browse files

Rotuv: added support for CDI_PROJ_RLL.

parent 9c598739
......@@ -33,22 +33,25 @@
static
void rot_uv_back(int gridID, double *us, double *vs)
{
long i, ilat, ilon, nlat, nlon;
long i, ilat, ilon;
double u, v;
double xval, yval;
double xpole, ypole, angle;
double *xvals, *yvals;
nlon = gridInqXsize(gridID);
nlat = gridInqYsize(gridID);
long nlon = gridInqXsize(gridID);
long nlat = gridInqYsize(gridID);
xpole = gridInqXpole(gridID);
ypole = gridInqYpole(gridID);
angle = gridInqAngle(gridID);
xvals = (double*) Malloc(nlon*sizeof(double));
yvals = (double*) Malloc(nlat*sizeof(double));
if ( gridInqType(gridID) == GRID_PROJECTION && gridInqProjType(gridID) == CDI_PROJ_RLL )
gridInqProjParamRLL(gridID, &xpole, &ypole, &angle);
else
{
xpole = gridInqXpole(gridID);
ypole = gridInqYpole(gridID);
angle = gridInqAngle(gridID);
}
double *xvals = (double*) Malloc(nlon*sizeof(double));
double *yvals = (double*) Malloc(nlat*sizeof(double));
gridInqXvals(gridID, xvals);
gridInqYvals(gridID, yvals);
......@@ -174,8 +177,9 @@ void *Rotuv(void *argument)
for ( varID = 0; varID < nvars; varID++ )
{
gridID = vlistInqVarGrid(vlistID1, varID);
if ( ! (gridInqType(gridID) == GRID_LONLAT && gridIsRotated(gridID)) )
cdoAbort("Only rotated lon/lat grids supported");
if ( ! ((gridInqType(gridID) == GRID_LONLAT && gridIsRotated(gridID)) ||
(gridInqType(gridID) == GRID_PROJECTION && gridInqProjType(gridID) == CDI_PROJ_RLL)) )
cdoAbort("Only rotated lon/lat grids supported!");
gridsize = gridInqSize(gridID);
nlevel = zaxisInqSize(vlistInqVarZaxis(vlistID1, varID));
......
......@@ -306,17 +306,13 @@ void grid_cell_center_to_bounds_Y2D(const char* yunitstr, long xsize, long ysize
}
}
void gridGenRotBounds(int gridID, int nx, int ny,
static
void gridGenRotBounds(double xpole, double ypole, double angle, int nx, int ny,
double *xbounds, double *ybounds, double *xbounds2D, double *ybounds2D)
{
double minlon, maxlon;
double minlat, maxlat;
double xpole = gridInqXpole(gridID);
double ypole = gridInqYpole(gridID);
double angle = gridInqAngle(gridID);
for ( int j = 0; j < ny; j++ )
{
if ( ybounds[0] > ybounds[1] )
......@@ -894,6 +890,7 @@ int gridToCurvilinear(int gridID1, int lbounds)
case GRID_LAEA:
case GRID_SINUSOIDAL:
{
double xpole = 0, ypole = 0, angle = 0;
double xscale = 1, yscale = 1;
double *xvals = NULL, *yvals = NULL;
double *xbounds = NULL, *ybounds = NULL;
......@@ -932,7 +929,7 @@ int gridToCurvilinear(int gridID1, int lbounds)
double *xvals2D = (double*) Malloc(gridsize*sizeof(double));
double *yvals2D = (double*) Malloc(gridsize*sizeof(double));
if ( gridtype == GRID_LCC )
if ( gridtype == GRID_LCC )
{
for ( int j = 0; j < ny; j++ )
for ( int i = 0; i < nx; i++ )
......@@ -967,8 +964,6 @@ int gridToCurvilinear(int gridID1, int lbounds)
if ( gridIsRotated(gridID1) || lrotated )
{
double xpole = 0, ypole = 0, angle = 0;
if ( lrotated )
{
gridInqProjParamRLL(gridID1, &xpole, &ypole, &angle);
......@@ -1103,9 +1098,9 @@ int gridToCurvilinear(int gridID1, int lbounds)
double *xbounds2D = (double*) Malloc(4*gridsize*sizeof(double));
double *ybounds2D = (double*) Malloc(4*gridsize*sizeof(double));
if ( gridIsRotated(gridID1) )
if ( gridIsRotated(gridID1) || lrotated )
{
gridGenRotBounds(gridID1, nx, ny, xbounds, ybounds, xbounds2D, ybounds2D);
gridGenRotBounds(xpole, ypole, angle, nx, ny, xbounds, ybounds, xbounds2D, ybounds2D);
}
else
{
......@@ -1308,6 +1303,7 @@ int gridToUnstructured(int gridID1, int lbounds)
case GRID_LONLAT:
case GRID_GAUSSIAN:
{
double xpole = 0, ypole = 0, angle = 0;
gridDefXname(gridID2, "lon");
gridDefYname(gridID2, "lat");
gridDefXlongname(gridID2, "longitude");
......@@ -1334,8 +1330,6 @@ int gridToUnstructured(int gridID1, int lbounds)
if ( gridIsRotated(gridID1) || lrotated )
{
double xpole = 0, ypole = 0, angle = 0;
if ( lrotated )
{
gridInqProjParamRLL(gridID1, &xpole, &ypole, &angle);
......@@ -1402,9 +1396,9 @@ int gridToUnstructured(int gridID1, int lbounds)
double *xbounds2D = (double*) Malloc(4*gridsize*sizeof(double));
double *ybounds2D = (double*) Malloc(4*gridsize*sizeof(double));
if ( gridIsRotated(gridID1) )
if ( gridIsRotated(gridID1) || lrotated )
{
gridGenRotBounds(gridID1, nx, ny, xbounds, ybounds, xbounds2D, ybounds2D);
gridGenRotBounds(xpole, ypole, angle, nx, ny, xbounds, ybounds, xbounds2D, ybounds2D);
}
else
{
......
Supports Markdown
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