PUBLIC INTERFACE ~ PUBLIC DATA ~ PUBLIC ROUTINES ~ NAMELIST ~ DIAGNOSTIC FIELDS ~ ERROR MESSAGES ~ REFERENCES ~ NOTES

Module horiz_interp_mod

Contact:  Bruce Wyman
Reviewers: 
Change History: WebCVS Log home


OVERVIEW

Performs spatial interpolation between grids.

This module can interpolate data from rectangular/tripolar grid to rectangular/tripolar grid. Three interpolation schems are used here. when the source grid is tripolar grid, inverse of square distance weighted scheme (spherical regrid) is used. When the source grid is rectangular longitude/latitude grid, any one of following three schemes can be applied: conservation scheme that conserves the area-weighed integral of the input field, bilinear interpolation that use the surround four source grid to interpolate onto the destination grid, spherical regrid that use thes inverse of square distance as weight. User can choose the interpolation method in the horiz_interp_init. The default method is conservative scheme. When the source grid is tripolar grid, only spherical regrid is allowed. When destination grid is tripolar, conservative scheme can not be used. There is an optional mask field for missing input data. An optional output mask field may be used in conjunction with the input mask to show where output data exists.


OTHER MODULES USED

      fms_mod
mpp_mod
constants_mod

PUBLIC INTERFACE

horiz_interp_init:
Allocates space and initializes a derived-type variable that contains pre-computed interpolation indices and weights.
horiz_interp:
Subroutine for performing the horizontal interpolation between two grids.
horiz_interp_end:
Deallocates memory used by "horiz_interp_type" variables. Must be called before reinitializing with horiz_interp_init.


PUBLIC DATA

None.


PUBLIC ROUTINES

  1. horiz_interp_init

     subroutine horiz_interp_init (Interp, lon_in, lat_in, lon_out, lat_out,  &
                                      verbose, interp_method, num_nbrs, max_dist, src_modulo)

    DESCRIPTION
    Allocates space and initializes a derived-type variable that contains pre-computed interpolation indices and weights for improved performance of multiple interpolations between the same grids. This routine does not need to be called if you are doing a single grid-to-grid interpolation.


    INPUT
    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 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 (when interp_method is "bilinear").
       [real, dimension(:),(:,:)]
    lat_out    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, 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]

    OUTPUT
    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)]

  2. horiz_interp

     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  )

    DESCRIPTION
    Subroutine for performing the horizontal interpolation between two grids. There are two forms of this interface. Form A requires first calling horiz_interp_init, while Form B requires no initialization.


    INPUT
    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(:),(:,:)]

    OUTPUT
    data_out    Output data on destination grid.
       [real, dimension(:,:),(:,:,:)]
    mask_out    Output mask that specifies whether data was computed.
       [real,optional, dimension(:,:),(:,:,:)]

  3. horiz_interp_end

    call horiz_interp_end ( Interp )
    DESCRIPTION
    Deallocates memory used by "horiz_interp_type" variables. Must be called before reinitializing with horiz_interp_init.


    INPUT/OUTPUT
    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]


PUBLIC TYPES

 type horiz_interp_type
   private
   real,    dimension(:,:), pointer   :: faci, facj  !weights for conservative scheme
   integer, dimension(:,:), pointer   :: ilon, jlat  !indices for conservative scheme
   real,    dimension(:,:), pointer   :: area_src    !area of the source grid
   real,    dimension(:,:), pointer   :: area_dst    !area of the destination grid
   real,    dimension(:,:,:), pointer :: wti,wtj     !weights for bilinear interpolation 
   integer, dimension(:,:,:), pointer :: i_lon, j_lat!indices for bilinear interpolation 
                                                     !and spherical regrid
   real,    dimension(:,:,:), pointer :: src_dist    !distance between destination grid and 
                                                     !neighbor source grid.
   logical, dimension(:,:), pointer   :: found_neighbors    !indicate whether destination grid 
                                                            !has some source grid around it.
   real                               :: max_src_dist
   integer                            :: nlon_src, nlat_src !size of source grid
   integer                            :: nlon_dst, nlat_dst !size of destination grid
   character(len=64)                  :: interp_method      !interpolation method.
                                                            !="conservative", conservative scheme
                                                            !="bilinear", bilinear interpolation
                                                            !="spherical", spherical regrid
 end type


DATA SETS

None.


ERROR MESSAGES

FATAL in horiz_interp
size of input array incorrect
The input data array does not match the size of the input grid edges specified. If you are using the initialization interface make sure you have the correct grid size.
FATAL in horiz_interp
size of output array incorrect
The output data array does not match the size of the input grid edges specified. If you are using the initialization interface make sure you have the correct grid size.


REFERENCES

None.


COMPILER SPECIFICS

None.


PRECOMPILER OPTIONS

None.


LOADER OPTIONS

None.


TEST PROGRAM

INTERNAL
       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


KNOWN BUGS

None.


NOTES

Has not been checked with grids that do not cover the sphere.

Has not been checked with the optional mask arguments.

If a latitude or longitude index cannot be found the tolerance used for making this determination may need to be increased. This can be done by increasing the value of module variable num_iters (default 4).


FUTURE PLANS

None.


top