Performs spatial interpolation between grids.
fms_mod
mpp_mod
constants_mod
horiz_interp_type_mod
horiz_interp_conserve_mod
horiz_interp_bilinear_mod
horiz_interp_spherical_mod
subroutine horiz_interp_init (Interp, lon_in, lat_in, lon_out, lat_out, & verbose, interp_method, num_nbrs, max_dist, src_modulo, & grid_at_center) |
lon_in | Longitude (in radians) for source data grid. when lon_in is 1D, it is the longitude
edges and its value are for adjacent grid boxes and must increase
in value. When lon_in is 2D, there are two cases: one is the longitude edges stored as
pairs for each grid box (when interp_method is "conservative"), the other is the longitude
of the center of each grid box. [real, dimension(:),(:,:)] |
lat_in | Latitude (in radians) for source data grid. when lat_in is 1D, it is the latitude
edges and its value are for adjacent grid boxes and must increase
in value. When lat_in is 2D, it is the longitude of the center of each grid box.
When interp_method is "conservative" or "bilinear", lon_in should be 1D. [real, dimension(:),(:,:)] |
lon_out | Longitude (in radians) for destination data grid. when lon_out is 1D, it is the longitude
edges and its value are for adjacent grid boxes and must increase
in value. When lon_out is 2D, there are two cases: one is the longitude edges stored as
pairs for each grid box (when interp_method is "conservative"), the other is the longitude
of the center of each grid box (when interp_method is "bilinear"). [real, dimension(:),(:,:)] |
lat_out | Latitude (in radians) for destination data grid. when lat_out is 1D, it is the latitude
edges and its value are for adjacent grid boxes and must increase
in value. When lat_out is 2D, there are two cases: one is the latitude edges stored as
pairs for each grid box (when interp_method is "conservative"), the other is the latitude
of the center of each grid box (when interp_method is "bilinear"). [real, dimension(:),(:,:)] |
verbose | Integer flag that controls the amount of printed output.
verbose = 0, no output; = 1, min,max,means; = 2, still more [integer, optional] |
interp_method | interpolation method, = "conservative", using conservation scheme,
= "bilinear", using bilinear interpolation, = "spherical",using spherical regrid.
when source grid is 1d, default value is "conservative"; when source grid is 2d,
default value is "spherical". [character(len=*),optional] |
src_modulo | Indicate the source data grid is cyclic or not. [logical, optional] |
grid_at_center | Indicate the data is on the center of grid box or the edge of grid box.
When true, the data is on the center of grid box. default vaule is false. [logical, optional] |
Interp | A derived-type variable containing indices and weights used for subsequent
interpolations. To reinitialize this variable for a different grid-to-grid
interpolation you must first use the "horiz_interp_end" interface. [type(horiz_interp_type)] |
subroutine horiz_interp ( Interp, data_in, data_out, verbose, & mask_in, mask_out, missing_value, missing_permit ) |
subroutine horiz_interp ( data_in, lon_in, lat_in, lon_out, lat_out, & data_out, verbose, mask_in, mask_out, & interp_method, missing_value, missing_permit, & num_nbrs, max_dist,src_modulo, grid_at_center ) |
Interp | Derived-type variable containing interpolation indices and weights.
Returned by a previous call to horiz_interp_init. [type(horiz_interp_type)] |
data_in | Input data on source grid. [real, dimension(:,:),(:,:,:)] |
verbose | flag for the amount of print output.
verbose = 0, no output; = 1, min,max,means; = 2, still more [integer,optional] |
mask_in | Input mask, must be the same size as the input data. The real value of
mask_in must be in the range (0.,1.). Set mask_in=0.0 for data points
that should not be used or have missing data. It is Not needed for
spherical regrid. [real,optional, dimension(:,:),(:,:,:)] |
missing_value | Use the missing_value to indicate missing data. [integer, optional] |
missing_permit | numbers of points allowed to miss for the bilinear interpolation. The value
should be between 0 and 3. [integer,optional] |
lon_in, lat_in | longitude and latitude (in radians) of source grid. More explanation can
be found in the documentation of horiz_interp_init. [real, dimension(:),(:,:)] |
lon_out, lat_out | longitude and latitude (in radians) of destination grid. More explanation can
be found in the documentation of horiz_interp_init. [real, dimension(:),(:,:)] |
data_out | Output data on destination grid. [real, dimension(:,:),(:,:,:)] |
mask_out | Output mask that specifies whether data was computed. [real,optional, dimension(:,:),(:,:,:)] |
call horiz_interp_end ( Interp )
Interp | A derived-type variable returned by previous call
to horiz_interp_init. The input variable must have
allocated arrays. The returned variable will contain
deallocated arrays. [horiz_interp_type] |
program test use horiz_interp_mod implicit none integer, parameter :: nxi=177, nyi=91, nxo=133, nyo=77 ! resolution real :: zi(nxi,nyi), zo(nxo,nyo) ! data real :: xi(nxi+1), yi(nyi+1), xo(nxo+1), yo(nyo+1) ! grid edges real :: pi, tpi, hpi, dx, dy constants hpi = acos(0.0) pi = hpi*2.0 tpi = hpi*4.0 grid setup: west to east, south to north dx = tpi/real(nxi); call setaxis (0.,dx,xi); xi(nxi+1) = xi(1)+tpi dx = tpi/real(nxo); call setaxis (0.,dx,xo); xo(nxo+1) = xo(1)+tpi dy = pi/real(nyi); call setaxis (-hpi,dy,yi); yi(nyi+1) = hpi dy = pi/real(nyo); call setaxis (-hpi,dy,yo); yo(nyo+1) = hpi random data on the input grid call random_number (zi) interpolate (flipping y-axis) call horiz_interp (zi(:,1:nyi:+1), xi, yi(1:nyi+1:+1), xo, yo(1:nyo+1:+1), zo, verbose=2) call horiz_interp (zi(:,nyi:1:-1), xi, yi(nyi+1:1:-1), xo, yo(1:nyo+1:+1), zo, verbose=2) call horiz_interp (zi(:,nyi:1:-1), xi, yi(nyi+1:1:-1), xo, yo(nyo+1:1:-1), zo, verbose=2) call horiz_interp (zi(:,1:nyi:+1), xi, yi(1:nyi+1:+1), xo, yo(nyo+1:1:-1), zo, verbose=2) contains set up a sequence of numbers subroutine setaxis (xo,dx,x) real, intent(in) :: xo, dx real, intent(out) :: x(:) integer :: i x(1) = xo do i=2,size(x(:)) x(i) = x(i-1)+dx enddo end subroutine setaxis end program test