INTERFACE:
subroutine areaini_point (io , jo , nlon_i , nlat_i , numlon_i, &
lon_i , lon_i_offset, lat_i , area_i , mask_i , &
nlon_o , nlat_o , numlon_o, lon_o , lat_o , &
area_o , fland_o , novr_i2o, iovr_i2o, jovr_i2o, &
wovr_i2o)
DESCRIPTION:
area averaging initialization
This subroutine is used for area-average mapping of a field from one
grid to another.
areaini_point - initializes indices and weights for area-averaging from
input grid to output grid
areamap_point - called by areaini: finds indices and weights
areaovr_point - called by areamap: finds if cells overlap and area of overlap
To map from one grid to another, must first call areaini to build
the indices and weights (iovr_i2o, jovr_i2o, wovr_i2o). Then must
call areaave to get new field on output grid.
Not all grid cells on the input grid will be used in the area-averaging
of a field to the output grid. Only input grid cells with [mask_i] = 1
contribute to output grid cell average. If [mask_i] = 0, input grid cell
does not contribute to output grid cell. This distinction is not usually
required for atm -> land mapping, because all cells on the atm grid have
data. But when going from land -> atm, only land grid cells have data.
Non-land grid cells on surface grid do not have data. So if output grid cell
overlaps with land and non-land cells (input grid), can only use land
grid cells when computing area-average.
o Input and output grids can be ANY resolution BUT:
a. Grid orientation -- Grids can be oriented south to north
(i.e. cell(lat+1) is north of cell(lat)) or from north to
south (i.e. cell(lat+1) is south of cell(lat)). Both grids must be
oriented from west to east, i.e., cell(lon+1) must be east of cell(lon)
b. Grid domain -- Grids do not have to be global. Both grids are defined
by their north, east, south, and west edges (edge_i and edge_o in
this order, i.e., edge_i(1) is north and edge_i(4) is west).
For partial grids, northern and southern edges are any latitude
between 90 (North Pole) and -90 (South Pole). Western and eastern
edges are any longitude between -180 and 180, with longitudes
west of Greenwich negative.
For global grids, northern and southern edges are 90 (North Pole)
and -90 (South Pole). The grids do not have to start at the
same longitude, i.e., one grid can start at Dateline and go east;
the other grid can start at Greenwich and go east. Longitudes for
the western edge of the cells must increase continuously and span
360 degrees. Examples
West edge East edge
---------------------------------------------------
Dateline : -180 to 180 (negative W of Greenwich)
Greenwich (centered): 0 - dx/2 to 360 - dx/2
c. Both grids can have variable number of longitude points for each
latitude strip. However, the western edge of the first point in each
latitude must be the same for all latitudes. Likewise, for the
eastern edge of the last point. That is, each latitude strip must span
the same longitudes, but the number of points to do this can be different
d. One grid can be a sub-set (i.e., smaller domain) than the other grid.
In this way, an atmospheric dataset for the entire globe can be
used in a simulation for a region 30N to 50N and 130W to 70W -- the
code will extract the appropriate data. The two grids do not have to
be the same resolution. Area-averaging will work for full => partial
grid but obviously will not work for partial => full grid.
o Field values fld_i on an input grid with dimensions nlon_i and nlat_i =>
field values fld_o on an output grid with dimensions nlon_o and nlat_o as
fld_o(io,jo) =
fld_i(i_ovr(io,jo, 1),j_ovr(io,jo, 1)) * w_ovr(io,jo, 1)
... + ... +
fld_i(i_ovr(io,jo,maxovr),j_ovr(io,jo,maxovr)) * w_ovr(io,jo,maxovr)
o Error checks:
Overlap weights of input cells sum to 1 for each output cell.
Global sum of dummy field is conserved for input => output area-average.
ARGUMENTS:
implicit none
integer , intent(in) :: io !output grid longitude index
integer , intent(in) :: jo !output grid latitude index
integer , intent(in) :: nlon_i !input grid: max number of longitude points
integer , intent(in) :: nlat_i !input grid: number of latitude points
integer , intent(in) :: numlon_i(nlat_i) !input grid: number lon points at each lat
real(r8), intent(in) :: lon_i(nlon_i+1,nlat_i) !input grid: longitude, west edge (degrees)
real(r8), intent(in) :: lon_i_offset(nlon_i+1,nlat_i) !input grid : cell lons, west edge (deg)
real(r8), intent(in) :: lat_i(nlat_i+1) !input grid: latitude, south edge (degrees)
real(r8), intent(in) :: area_i(nlon_i,nlat_i) !input grid: cell area
real(r8), intent(in) :: mask_i(nlon_i,nlat_i) !input grid: mask (0, 1)
integer , intent(in) :: nlon_o !output grid: max number of longitude points
integer , intent(in) :: nlat_o !output grid: number of latitude points
integer , intent(in) :: numlon_o(nlat_o) !output grid: number lon points at each lat
real(r8), intent(in) :: lon_o(nlon_o+1,nlat_o) !output grid: longitude, west edge (degrees)
real(r8), intent(in) :: lat_o(nlat_o+1) !output grid: latitude, south edge (degrees)
real(r8), intent(in) :: area_o !output grid: cell area
real(r8), intent(in) :: fland_o !output grid: fraction that is land
integer , intent(out) :: novr_i2o !number of overlapping input cells
integer , intent(out) :: iovr_i2o(maxovr) !lon index of overlap input cell
integer , intent(out) :: jovr_i2o(maxovr) !lat index of overlap input cell
real(r8), intent(out) :: wovr_i2o(maxovr) !weight of overlap input cell
REVISION HISTORY:
Created by Gordon Bonan