3.2.2. Grids
[3.2. Reading Geometry]

The XMDF library supports many different types of grids. Grids can be two or three dimensional. Two-dimensional grids may be Cartesian or curvilinear. Three-dimensional grids may be Cartesian, curvilinear, or extruded 2D Cartesian or curvilinear grids. Extruded grids may be sigma-stretch, curvilinear, curvilinear at corners, or curvilinear at mid-sides. Grid properties include the origin, orientation, dip, and bearing. The way the grid geometry is specified depends upon the type of grid. For a Cartesian grid, the gridline locations along each axis (I, J, and K if applicable) must be specified. These locations are defined as the distance along the axis from the grid origin. For a curvilinear grid, the locations of each grid corner must be defined. The locations are given in world coordinates.

To read a grid from an XMDF file:

Open the file and the grid group using xfOpenFile and xfOpenGroup.

Make sure that the grid type and properties are valid for the particular model. If the model requires a 3D curvilinear grid, you need to abort if the grid is any other type. Functions that may be used include xfGetGridType, xfGetExtrusionType, xfGetNumberOfDimensions.

Get the number of cells in each direction using xfGetNumberCellsInI, xfGetNumberCellsInJ, and xfGetNumberCellsInK.

Allocate the arrays for the grid geometry definition. Remember that the size of the arrays depends not only upon the number of cells in the grid but also the type of grid.

Read the grid geometry using xfGetGridCoordsI, xfGetGridCoordsJ, and xfGetGridCoordsK.

Close the file and grid group using xfCloseFile and xfCloseGroup.

Convert grid information into native model data definitions if necessary.

Example

The example below has a function to read data for a model that requires a 3D Cartesian grid for input. The model uses a bearing but dip and roll values are ignored. The filename and path to the grid is given and the data is stored in variables accessible outside the function (file globals in C/C++, common blocks in FORTRAN). The function follows the XMDF convention of using negative values to indicate errors.

C/C++
int       fg_nCellsI, fg_nCellsJ, fg_nCellsK;
double    fg_Origin[3], fg_Bearing;
double   *fg_CoordI, *fg_CoordJ, *fg_CoordK;
 
int ReadGrid (const char *a_Filename, const char *a_GridPath)
{
  xid      xFileId, xGroupId;
  int      nGridType = 0, nDims = 0;
  int      nValsI = 0, nValsJ = 0;
  int      nValsK = 0;
  xbool    bDefined = XFALSE;
  int      i = 0, error = 1;
 
  // open the file and group
  error = xfOpenFile(a_Filename, &xFileId, TRUE);
  if (error < 0) {
    printf("Unable to open file.");
    return -1;
  }
 
  error = xfOpenGroup(xFileId, a_GridPath, &xGroupId);
  if (error < 0) {
    printf("Unable to open the group.");
    xfCloseFile(xFileId);
    return -1;
  }
 
  // Grid type
  error = xfGetGridType(a_Id, &nGridType);
  if (error < 0) {
    return error;
  }
  if (nGridType != GRID_TYPE_CARTESIAN) {
    printf("Unsupported grid type. Must be a Cartesian Grid");
    return -1;
  }
  // Number of dimensions
  error = xfGetNumberOfDimensions(a_Id, &nDims);
  if (error < 0) {
    return error;
  }
  if (nDims != 3) {
    printf("Error:  The grid must be a three-dimensional grid.");
    return -1;
  }
 
  // Origin
  error = xfOriginDefined(a_Id, &bDefined);
  if (error < 0) {
    return error;
  }
  if (bDefined) {
    error = xfGetOrigin(a_Id, &fg_Origin[0], &fg_Origin[1], &fg_Origin[2]);
    if (error < 0) {
      return error;
    }
  }
  
  // Bearing
  error = xfBearingDefined(a_Id, &bDefined);
  if (error < 0) {
    return error;
  }
  if (bDefined) {
    error = xfGetBearing(a_Id, &fg_Bearing);
    if (error < 0) {
      return error;
    }
  }
 
  // number of cells in each direction
  error = xfGetNumberCellsInI(a_Id, &nCellsI);
  if (error >= 0) {
    error = xfGetNumberCellsInJ(a_Id, &nCellsJ);
    if (error >= 0 && nDims == 3) { 
      error = xfGetNumberCellsInK(a_Id, &nCellsK);
    }
  }
  if (error < 0) {
    return error;
  }
 
  nValsI = fg_nCellsI;
  nValsJ = fg_nCellsJ;
  nValsK = fg_nCellsK;
 
  fg_CoordI = new double[nValsI];
  fg_CoordJ = new double[nValsJ];
  fg_CoordK = new double[nValsK];
 
  error = xfGetGridCoordsI(a_Id, nValsI, fg_CoordI);
  if (error >= 0) {
    error = xfGetGridCoordsJ(a_Id, nValsJ, fg_CoordJ);
    if (error >= 0 && nDims == 3) {
      error = xfGetGridCoordsK(a_Id, nValsK, fg_CoordK);
    }
  }
  if (error < 0) {
    printf("Error reading coordinates.\n");
    return -1;
  }
 
  // if we got here, everything read in fine
  return 1;
} // ReadGrid

FORTRAN
INTEGER      fg_nCellsI, fg_nCellsJ, fg_nCellsK
REAL(DOUBLE), DIMENSION(3) :: fg_Origin
REAL(DOUBLE), ALLOCATABLE :: fg_CoordI(:), fg_CoordJ(:), fg_CoordK(:)
REAL(DOUBLE)              :: fg_Bearing
 
SUBROUTINE  READ_GRID(a_Filename, a_GridPath, error)
CHARACTER(LEN=*), INTENT(IN) :: a_Filename, a_GridPath 
INTEGER, INTENT(OUT)       :: error
INTEGER(XID)   xFileId, xGroupId
INTEGER nGridType, nDims;
INTEGER nValsI, nValsJ, nValsK
INTEGER bDefined
 
call XF_OPEN_FILE(a_Filename, XTRUE, xFileId, error)
if (error < 0) then
  return
endif
 
call XF_OPEN_GROUP(xFileId, a_GridPath, xGroupId, error)
if (error < 0) then
  XF_CLOSE_FILE(xFileId, error)
  return
endif
 
  ! Grid type
call XF_GET_GRID_TYPE(a_Id, nGridType, error)
if (error < 0) then
  return
endif
if (nGridType .EQ. GRID_TYPE_CARTESIAN) then
  Write(*,*) 'Unsupported grid type.  Must be a Cartesian Grid'
  return
endif
 
  ! Number of dimensions
call XF_GET_NUMBER_OF_DIMENSIONS(a_Id, nDims, error)
if (error .LT. 0) then
  return
endif
if (nDims .NE. 3) then
  WRITE(*,*) 'The grid must be a three-dimensional grid'
  error = -1
  return
endif
 
  ! Origin
call XF_ORIGIN_DEFINED(a_Id, bDefined, error)
if (error < 0) then
  return
endif
if (bDefined /= 0) then
  call XF_GET_ORIGIN(a_Id, fg_Origin(1), fg_Origin(2), fg_Origin(3), error)
  if (error < 0) then
    return
  endif
endif
  
  ! Bearing
call XF_BEARING_DEFINED(a_Id, bDefined, error)
if (error < 0) then
  return
endif
if (bDefined /= 0) then
  call XF_GET_BEARING(a_Id, fg_Bearing, error)
  if (error < 0) then
    return
  endif
endif
 
  ! number of cells in each direction
call XF_GET_NUMBER_CELLS_IN_I(a_Id, fg_nCellsI, error)
if (error >= 0) then
  call XF_GET_NUMBER_CELLS_IN_J(a_Id, fg_nCellsJ, error)
  if (error > 0) then 
    call XF_GET_NUMBER_CELLS_IN_K(a_Id, fg_nCellsK, error)
  endif
endif
if (error < 0) then
  return
endif
 
nValsI = fg_nCellsI
nValsJ = fg_nCellsJ
nValsK = fg_nCellsK
 
ALLOCATE(fg_CoordI(nValsI))
ALLOCATE(fg_CoordJ(nValsJ))
ALLOCATE(fg_CoordK(nValsK))
 
call XF_GET_GRID_COORDS_I(a_Id, nValsI, fg_CoordI, error)
if (error >= 0) then
  call XF_GET_GRID_COORDS_J(a_Id, nValsJ, fg_CoordJ, error)
  if ((error > 0) .AND. (nDims == 3)) then
    call XF_GET_GRID_COORDS_K(a_Id, nValsK, fg_CoordK, error)
  endif
endif
if (error < 0) then
  WRITE(*,*) 'Error reading coordinates'
  error = -1
  return
endif
 
! Return successful if we got here
error = 1
return
END SUBROUTINE TG_READ_GRID