FMS 2025.01-dev
Flexible Modeling System
Loading...
Searching...
No Matches
diag_grid_mod

diag_grid_mod is a set of procedures to work with the model's global grid to allow regional output. More...

Data Types

type  diag_global_grid_type
 Private type to hold the model's global grid data, and other grid information for use in this module. More...
 
type  point
 Private type to hold the corresponding (x,y,z) location for a (lat,lon) location. More...
 

Functions/Subroutines

pure elemental real function deg2rad (angle)
 Convert an angle in degrees to radians.
 
subroutine, public diag_grid_end ()
 Unallocate the diag_global_grid variable.
 
subroutine, public diag_grid_init (domain, glo_lat, glo_lon, aglo_lat, aglo_lon)
 Send the global grid to the diag_manager_mod for regional output.
 
pure elemental real function distancesqrd (pt1, pt2)
 Find the distance between two points in the Cartesian coordinate space.
 
pure integer function, dimension(2) find_equator_index_agrid (lat, lon)
 Return the closest index (i,j) to the given (lat,lon) point.
 
subroutine find_nearest_agrid_index (lat, lon, mini, minj, minimum_distance)
 Find the i,j indices and distance of the a-grid point nearest to the inputted lat,lon point.
 
pure integer function, dimension(2) find_pole_index_agrid (lat, lon)
 Return the closest index (i,j) to the given (lat,lon) point.
 
pure elemental real function gcirdistance (lat1, lon1, lat2, lon2)
 Find the distance, along the geodesic, between two points.
 
subroutine, public get_local_indexes (latstart, latend, lonstart, lonend, istart, iend, jstart, jend)
 Find the local start and local end indexes on the local PE for regional output.
 
subroutine, public get_local_indexes2 (lat, lon, iindex, jindex)
 Find the indices of the nearest grid point of the a-grid to the specified (lon,lat) location on the local PE. if desired point not within domain of local PE, return (0,0) as the indices.
 
pure elemental type(point) function latlon2xyz (lat, lon)
 Return the (x,y,z) position of a given (lat,lon) point.
 
pure elemental real function rad2deg (angle)
 

Variables

type(diag_global_grid_typediag_global_grid
 Variable to hold the global grid data.
 
logical diag_grid_initialized = .FALSE.
 Indicates if the diag_grid_mod has been initialized.
 

Detailed Description

diag_grid_mod is a set of procedures to work with the model's global grid to allow regional output.

Author
Seth Underwood seth..nosp@m.unde.nosp@m.rwood.nosp@m.@noa.nosp@m.a.gov

diag_grid_mod contains useful utilities for dealing with, mostly, regional output for grids other than the standard lat/lon grid. This module contains three public procedures diag_grid_init, which is shared globably in the diag_manager_mod, diag_grid_end which will free up memory used during the register field calls, and get_local_indexes. The send_global_grid procedure is called by the model that creates the global grid. send_global_grid needs to be called before any fields are registered that will output only regions. get_local_indexes is to be called by the diag_manager_mod to discover the global indexes defining a subregion on the tile.


Data Type Documentation

◆ diag_grid_mod::diag_global_grid_type

type diag_grid_mod::diag_global_grid_type

Private type to hold the model's global grid data, and other grid information for use in this module.

Definition at line 73 of file diag_grid.F90.

Collaboration diagram for diag_global_grid_type:
[legend]

Public Attributes

integer adimi
 The dimension of the global a-grid in the 'i' / longitudal direction. Again, the expected dimension for diag_grid_mod is isc-1:iec+1.
 
integer adimj
 The dimension of the global a-grid in the 'j' / latitudal direction. Again, the expected dimension for diag_grid_mod is jsc-1:jec+1.
 
real, dimension(:,:), allocatable aglo_lat
 The latitude values on the global a-grid. Here we expect isc-1:iec+1 and jsc=1:jec+1 to be passed in.
 
real, dimension(:,:), allocatable aglo_lon
 The longitude values on the global a-grid. Here we expec isc-1:iec+j and jsc-1:jec+1 to be passed in.
 
integer dimi
 The dimension of the global grid in the 'i' / longitudal direction.
 
integer dimj
 The dimension of the global grid in the 'j' / latitudal direction.
 
real, dimension(:,:), allocatable glo_lat
 The latitude values on the global grid.
 
real, dimension(:,:), allocatable glo_lon
 The longitude values on the global grid.
 
character(len=128) grid_type
 The global grid type.
 
integer myxbegin
 The starting index of the compute domain on the current PE.
 
integer myybegin
 The starting index of the compute domain on the current PE.
 
integer ntiles
 The number of tiles.
 
integer peend
 The ending PE number for the current tile.
 
integer pestart
 The starting PE number for the current tile.
 
integer tile_number
 The tile the glo_lat and glo_lon define.
 

Member Data Documentation

◆ adimi

integer adimi

The dimension of the global a-grid in the 'i' / longitudal direction. Again, the expected dimension for diag_grid_mod is isc-1:iec+1.

Definition at line 86 of file diag_grid.F90.

◆ adimj

integer adimj

The dimension of the global a-grid in the 'j' / latitudal direction. Again, the expected dimension for diag_grid_mod is jsc-1:jec+1.

Definition at line 88 of file diag_grid.F90.

◆ aglo_lat

real, dimension(:,:), allocatable aglo_lat

The latitude values on the global a-grid. Here we expect isc-1:iec+1 and jsc=1:jec+1 to be passed in.

Definition at line 76 of file diag_grid.F90.

◆ aglo_lon

real, dimension(:,:), allocatable aglo_lon

The longitude values on the global a-grid. Here we expec isc-1:iec+j and jsc-1:jec+1 to be passed in.

Definition at line 79 of file diag_grid.F90.

◆ dimi

integer dimi

The dimension of the global grid in the 'i' / longitudal direction.

Definition at line 84 of file diag_grid.F90.

◆ dimj

integer dimj

The dimension of the global grid in the 'j' / latitudal direction.

Definition at line 85 of file diag_grid.F90.

◆ glo_lat

real, dimension(:,:), allocatable glo_lat

The latitude values on the global grid.

Definition at line 74 of file diag_grid.F90.

◆ glo_lon

real, dimension(:,:), allocatable glo_lon

The longitude values on the global grid.

Definition at line 75 of file diag_grid.F90.

◆ grid_type

character(len=128) grid_type

The global grid type.

Definition at line 94 of file diag_grid.F90.

◆ myxbegin

integer myxbegin

The starting index of the compute domain on the current PE.

Definition at line 82 of file diag_grid.F90.

◆ myybegin

integer myybegin

The starting index of the compute domain on the current PE.

Definition at line 83 of file diag_grid.F90.

◆ ntiles

integer ntiles

The number of tiles.

Definition at line 91 of file diag_grid.F90.

◆ peend

integer peend

The ending PE number for the current tile.

Definition at line 93 of file diag_grid.F90.

◆ pestart

integer pestart

The starting PE number for the current tile.

Definition at line 92 of file diag_grid.F90.

◆ tile_number

integer tile_number

The tile the glo_lat and glo_lon define.

Definition at line 90 of file diag_grid.F90.

◆ diag_grid_mod::point

type diag_grid_mod::point

Private type to hold the corresponding (x,y,z) location for a (lat,lon) location.

Definition at line 100 of file diag_grid.F90.

Collaboration diagram for point:
[legend]

Public Attributes

real x
 The x value of the (x,y,z) coordinates.
 
real y
 The y value of the (x,y,z) coordinates.
 
real z
 The z value of the (x,y,z) coordinates.
 

Member Data Documentation

◆ x

real x

The x value of the (x,y,z) coordinates.

Definition at line 101 of file diag_grid.F90.

◆ y

real y

The y value of the (x,y,z) coordinates.

Definition at line 102 of file diag_grid.F90.

◆ z

real z

The z value of the (x,y,z) coordinates.

Definition at line 103 of file diag_grid.F90.

Function/Subroutine Documentation

◆ deg2rad()

pure elemental real function deg2rad ( real, intent(in)  angle)
private

Convert an angle in degrees to radians.

Given a scalar, or an array of angles in degrees this function will return a scalar or array (of the same dimension) of angles in radians.

Returns
Scalar or array (depending on the size of angle) of angles in radians.
Parameters
[in]angleScalar or array of angles in degrees.

Definition at line 696 of file diag_grid.F90.

◆ diag_grid_end()

subroutine, public diag_grid_end

Unallocate the diag_global_grid variable.

The diag_global_grid variable is only needed during the register field calls, and then only if there are fields requestion regional output. Once all the register fields calls are complete (before the first send_data call this procedure can be called to free up memory.

Definition at line 349 of file diag_grid.F90.

◆ diag_grid_init()

subroutine, public diag_grid_init ( type(domain2d), intent(in)  domain,
class(*), dimension(:,:), intent(in)  glo_lat,
class(*), dimension(:,:), intent(in)  glo_lon,
class(*), dimension(:,:), intent(in)  aglo_lat,
class(*), dimension(:,:), intent(in)  aglo_lon 
)

Send the global grid to the diag_manager_mod for regional output.

In order for the diag_manager to do regional output for grids other than the standard lat/lon grid, the diag_manager_mod needs to know the the latitude and longitude values for the entire global grid. This procedure is the mechanism the models will use to share their grid with the diagnostic manager.

This procedure needs to be called after the grid is created, and before the first call to register the fields.

Parameters
[in]domainThe domain to which the grid data corresponds.
[in]glo_latThe latitude information for the grid tile.
[in]glo_lonThe longitude information for the grid tile.
[in]aglo_latThe latitude information for the a-grid tile.
[in]aglo_lonThe longitude information for the a-grid tile.

Definition at line 131 of file diag_grid.F90.

◆ distancesqrd()

pure elemental real function distancesqrd ( type(point), intent(in)  pt1,
type(point), intent(in)  pt2 
)
private

Find the distance between two points in the Cartesian coordinate space.

distanceSqrd will find the distance squared between two points in the xyz coordinate space. pt1 and pt2 can either be both scalars, both arrays of the same size, or one a scalar and one an array. The return value will be a scalar or array of the same size as the input array.

Returns
The return value will be a scalar or array of the same size as the input array.

Definition at line 1054 of file diag_grid.F90.

◆ find_equator_index_agrid()

pure integer function, dimension(2) find_equator_index_agrid ( real, intent(in)  lat,
real, intent(in)  lon 
)
private

Return the closest index (i,j) to the given (lat,lon) point.

This function searches a equator grid tile looking for the grid point closest to the give (lat, lon) location, and returns the i and j indexes of the point.

Returns
The (i, j) location of the closest grid to the given (lat, lon) location.
The (i, j) location of the closest grid to the given (lat, lon) location.
Parameters
[in]latLatitude location
[in]lonLongitude location

Definition at line 873 of file diag_grid.F90.

◆ find_nearest_agrid_index()

subroutine find_nearest_agrid_index ( real, intent(in)  lat,
real, intent(in)  lon,
integer, intent(out)  mini,
integer, intent(out)  minj,
real, intent(out)  minimum_distance 
)
private

Find the i,j indices and distance of the a-grid point nearest to the inputted lat,lon point.

Definition at line 1084 of file diag_grid.F90.

◆ find_pole_index_agrid()

pure integer function, dimension(2) find_pole_index_agrid ( real, intent(in)  lat,
real, intent(in)  lon 
)
private

Return the closest index (i,j) to the given (lat,lon) point.

This function searches a pole a-grid tile looking for the grid point closest to the give (lat, lon) location, and returns the i and j indexes of the point.

Returns
The (i, j) location of the closest grid to the given (lat, lon) location.
The (i, j) location of the closest grid to the given (lat, lon) location.
Parameters
[in]latLatitude location
[in]lonLongitude location

Definition at line 709 of file diag_grid.F90.

◆ gcirdistance()

pure elemental real function gcirdistance ( real, intent(in)  lat1,
real, intent(in)  lon1,
real, intent(in)  lat2,
real, intent(in)  lon2 
)
private

Find the distance, along the geodesic, between two points.

gCirDistance will find the distance, along the geodesic, between two points defined by the (lat,lon) position of each point.

Returns
real

Definition at line 1067 of file diag_grid.F90.

◆ get_local_indexes()

subroutine, public get_local_indexes ( real, intent(in)  latstart,
real, intent(in)  latend,
real, intent(in)  lonstart,
real, intent(in)  lonend,
integer, intent(out)  istart,
integer, intent(out)  iend,
integer, intent(out)  jstart,
integer, intent(out)  jend 
)

Find the local start and local end indexes on the local PE for regional output.

Given a defined region, find the local indexes on the local PE surrounding the region.

Parameters
[in]latstartlat start angles
[in]lonstartlon start angles
[in]latendlat end angles
[in]lonendlon end angles
[out]istarti start indexes
[out]jstartj start indexes
[out]iendi end indexes
[out]jendj end indexes

Definition at line 390 of file diag_grid.F90.

◆ get_local_indexes2()

subroutine, public get_local_indexes2 ( real, intent(in)  lat,
real, intent(in)  lon,
integer, intent(out)  iindex,
integer, intent(out)  jindex 
)

Find the indices of the nearest grid point of the a-grid to the specified (lon,lat) location on the local PE. if desired point not within domain of local PE, return (0,0) as the indices.

Parameters
[in]latlat location
[in]lonlon location
[out]iindexi indexes
[out]jindexj indexes

Definition at line 640 of file diag_grid.F90.

◆ latlon2xyz()

pure elemental type(point) function latlon2xyz ( real, intent(in)  lat,
real, intent(in)  lon 
)
private

Return the (x,y,z) position of a given (lat,lon) point.

Given a specific (lat, lon) point on the Earth, return the corresponding (x,y,z) location. The return of latlon2xyz will be either a scalar or an array of the same size as lat and lon.

Returns
The return of latlon2xyz will be either a scalar or an array of the same size as lat and lon.
Parameters
[in]latThe latitude of the (x,y,z) location to find. lat can be either a scalar or array. lat must be of the same rank / size as lon. This function assumes lat is in the range [-90,90].
[in]lonThe longitude of the (x,y,z) location to find. lon can be either a scalar or array. lon must be of the same rank / size as lat. This function assumes lon is in the range [0,360].

Definition at line 1018 of file diag_grid.F90.

◆ rad2deg()

pure elemental real function rad2deg ( real, intent(in)  angle)
private
Parameters
[in]angleScalar or array of angles in radians.

Definition at line 683 of file diag_grid.F90.

Variable Documentation

◆ diag_global_grid

type(diag_global_grid_type) diag_global_grid

Variable to hold the global grid data.

Definition at line 109 of file diag_grid.F90.

◆ diag_grid_initialized

logical diag_grid_initialized = .FALSE.

Indicates if the diag_grid_mod has been initialized.

Definition at line 111 of file diag_grid.F90.