User's Guide

SURFACE

Grid elevation-dependent observations from scattered stations

Function:

Surface creates a uniformly spaced surface grid from values of an elevation-dependent variable recorded at a set of observing stations which are not uniformly distributed. The contribution of stations near a given grid point are weighted by an exponential function of the distance from the point to the station, subject to a cutoff radius. The value of the variable at each contributing observing station is adjusted to the elevation of the grid point using the slope and intercept of a regression line separately computed from, typically, the full set of stations.

This program does not apply spatial and temporal smoothing to observation and station elevation data -- see User Note 1.

The surfacing algorithm is taken from P. E. Thornton, S. W. Running, and M. A. White, "Generating surfaces of daily meteorological variables over large regions of complex terrain", Journal of Hydrology 190 (1997), pp. 214-251.

Parameters:

Subcommand -GENERAL:
Surface an arbitrary variable. The default values of the arguments are the same as for subcommand TEMP.

INELEV
Image containing elevation for each output grid point. Values should be in meters. This single-band image defines the projection and spatial extent of the output image; standard LAS window specifiers are supported. The INELEV image could also specify some other surface property on which the observed values depend -- see Usernote 3.

INOBS
File containing location of each observing station and value of variable to be surfaced. There must be one record for each station, containing four whitespace-separated fields giving the latitude and longitude, in decimal degrees, the elevation in meters, and the value of the observed variable as a real or integer value. South latitudes and west longitudes must be indicated by negative values.

A constant value for the slope of the regression line between observed values and elevation may be specified in the first record of this file, in the format "slope=XXX.XXXX", where the number of places on each side of the decimal point is arbitrary; this value is ignored if argument SLOPE is specified.

OUT
Name of output image. Will contain gridded, elevation adjusted values of observed varialbe.

ODTYPE("R*4")
Data type of output image.

ALPHA(3.0)
Shape parameter of Gaussian weighting function.

STACNT(30)
Target maximum number of stations for each grid point. The cutoff distance for the Gaussian weighting function will be adjusted for each region of the image to include approximately this number of stations. In practice, due to differenct patterns of station distribution, the actual number of stations used for a given grid point may be somewhat less.

NITER(3)
Number of iterations for refining maximum distance. An initial trial value of the cutoff distance is computed from the mean density of observation stations over the entire image window. Since the station density may vary significantly over the image, the cutoff distance is refined iteratively for each image region.

SLOPE(--)
Constant value for slope of regression line between variable and elevation. This overrides a slope value specified in the INOBS file. If SLOPE is defaulted and no value is specified in the INOBS file, a weighted regression is computed for each output pixel.

REGRESS("LINEAR")
Type of regression of variable against elevation. Options are


    NONE:	Do not adjust value of variable for elevation.
    LINEAR:	Use linear relation between variable and elevation.
    NORMDIFF:	Assume linear relation between elevation and the
		normalized difference between the value V of the
		variable at a given grid point and its value Vo at the
		observing station; ND = (V - Vo) / (V + Vo).

MAXND(0.6)
Maximum absolute value of normalized difference. Normalized differences exceeding MAXND will be replaced by MAXND, with the same sign. This applies both to normalized difference values between pairs of stations used for computing the regression coefficients, and for the values used to compute output pixels. The default value, 0.6, corresponds to a maximum ratio between values of 4 to 1.

MINFRACT(0.0)
Minimum fraction of stations above threshold for non-fill output. If MINFRACT is greater than zero, the program determines the weighted fraction of observing stations contributing to each output pixel for which the observed variable has a value greater than or equal to THRESH. If this fraction is less than MINFRACT, the output pixel's value will be the value specified by FILL.

THRESH(0.0)
Threshold for using variable value from station. If the observed value of the variable at the station is below THRESH, the station will not be used for computing the weighted value of the output pixel.

FILL(0.0)
Fill value for output pixels failing MINFRACT test.

PRINT(TERM)
Destination(s) for processing report. This summarizes the average, maximum, and minimum sum of the processing weights for each grid cell, the adjusted cutoff distance from a point to a station, and the values of slope and intercept used. This listing may be sent to one or more destinations as specifed by the value(s) entered for PRINT. Options are


    TERM:	Listing is displayed on the user's terminal.
    LP:		Listing is sent to the  printer defined by $PRINTER.
    filename:	Listing is stored in a file with the specified filename;
		the file extension ".prt" is appended to the filename
		specified by the user.
Subcommand -TEMP:
Surface temperature values. The default values of the arguments specify a linear regression between the observation data and elevation, using the weighting parameters and number of contributing stations which were found to be optimum for the Pacific Northwest in the paper by Thornton, et al.

INELEV
Image containing elevation for each output grid point. Values should be in meters. This single-band image defines the projection and spatial extent of the output image; standard LAS window specifiers are supported. The INELEV image could also specify some other surface property on which the observed values depend -- see Usernote 3.

INOBS
File containing location of each observing station and value of variable to be surfaced. There must be one record for each station, containing four whitespace-separated fields giving the latitude and longitude, in decimal degrees, the elevation in meters, and the value of the observed variable as a real or integer value. South latitudes and west longitudes must be indicated by negative values.

A constant value for the slope of the regression line between observed values and elevation may be specified in the first record of this file, in the format "slope=XXX.XXXX", where the number of places on each side of the decimal point is arbitrary; this value is ignored if argument SLOPE is specified.

OUT
Name of output image. Will contain gridded, elevation adjusted values of observed varialbe.

ODTYPE("BYTE")
Data type of output image.

ALPHA(3.0)
Shape parameter of Gaussian weighting function.

STACNT(30)
Target maximum number of stations for each grid point. The cutoff distance for the Gaussian weighting function will be adjusted for each region of the image to include approximately this number of stations. In practice, due to differenct patterns of station distribution, the actual number of stations used for a given grid point may be somewhat less.

NITER(3)
Number of iterations for refining maximum distance. An initial trial value of the cutoff distance is computed from the mean density of observation stations over the entire image window. Since the station density may vary significantly over the image, the cutoff distance is refined iteratively for each image region.

SLOPE(--)
Constant value for slope of regression line between variable and elevation. This overrides a slope value specified in the INOBS file. If SLOPE is defaulted and no value is specified in the INOBS file, a weighted regression is computed for each output pixel.

REGRESS("LINEAR")
Type of regression of variable against elevation. Options are


    NONE:	Do not adjust value of variable for elevation.
    LINEAR:	Use linear relation between variable and elevation.
    NORMDIFF:	Assume linear relation between elevation and the
		normalized difference between the value V of the
		variable at a given grid point and its value Vo at the
		observing station; ND = (V - Vo) / (V + Vo).

MAXND(0.6)
Maximum absolute value of normalized difference. Normalized differences exceeding MAXND will be replaced by MAXND, with the same sign. This applies both to normalized difference values between pairs of stations used for computing the regression coefficients, and for the values used to compute output pixels. The default value, 0.6, corresponds to a maximum ratio between values of 4 to 1.

MINFRACT(0.0)
Minimum fraction of stations above threshold for non-fill output. If MINFRACT is greater than zero, the program determines the weighted fraction of observing stations contributing to each output pixel for which the observed variable has a value greater than or equal to THRESH. If this fraction is less than MINFRACT, the output pixel's value will be the value specified by FILL.

THRESH(0.0)
Threshold for using variable value from station. If the observed value of the variable at the station is below THRESH, the station will not be used for computing the weighted value of the output pixel.

FILL(0.0)
Fill value for output pixels failing MINFRACT test.

PRINT(TERM)
Destination(s) for processing report. This summarizes the average, maximum, and minimum sum of the processing weights for each grid cell, the adjusted cutoff distance from a point to a station, and the values of slope and intercept used. This listing may be sent to one or more destinations as specifed by the value(s) entered for PRINT. Options are


    TERM:	Listing is displayed on the user's terminal.
    LP:		Listing is sent to the  printer defined by $PRINTER.
    filename:	Listing is stored in a file with the specified filename;
		the file extension ".prt" is appended to the filename
		specified by the user.
Subcommand -PRECIP:
Surface precipitation values. The default values of the arguments specify a regression of the normed difference between observations against station elevation, using the weighting parameters and number of contributing stations which were found to be optimum for the Pacific Northwest in the paper by Thornton, et al. The program does not perform the spatial and temporal smoothing recommended in that paper -- see Usernote 1.

INELEV
Image containing elevation for each output grid point. Values should be in meters. This single-band image defines the projection and spatial extent of the output image; standard LAS window specifiers are supported. The INELEV image could also specify some other surface property on which the observed values depend -- see Usernote 3.

INOBS
File containing location of each observing station and value of variable to be surfaced. There must be one record for each station, containing four whitespace-separated fields giving the latitude and longitude, in decimal degrees, the elevation in meters, and the value of the observed variable as a real or integer value. South latitudes and west longitudes must be indicated by negative values.

A constant value for the slope of the regression line between observed values and elevation may be specified in the first record of this file, in the format "slope=XXX.XXXX", where the number of places on each side of the decimal point is arbitrary; this value is ignored if argument SLOPE is specified.

OUT
Name of output image. Will contain gridded, elevation adjusted values of observed varialbe.

ODTYPE("R*4")
Data type of output image.

ALPHA(6.25)
Shape parameter of Gaussian weighting function.

STACNT(20)
Target maximum number of stations for each grid point. The cutoff distance for the Gaussian weighting function will be adjusted for each region of the image to include approximately this number of stations. In practice, due to differenct patterns of station distribution, the actual number of stations used for a given grid point may be somewhat less.

NITER(3)
Number of iterations for refining maximum distance. An initial trial value of the cutoff distance is computed from the mean density of observation stations over the entire image window. Since the station density may vary significantly over the image, the cutoff distance is refined iteratively for each image region.

SLOPE(--)
Constant value for slope of regression line between variable and elevation. This overrides a slope value specified in the INOBS file. If SLOPE is defaulted and no value is specified in the INOBS file, a weighted regression is computed for each output pixel.

REGRESS("NORMDIFF")
Type of regression of variable against elevation. Options are


    NONE:	Do not adjust value of variable for elevation.
    LINEAR:	Use linear relation between variable and elevation.
    NORMDIFF:	Assume linear relation between elevation and the
		normalized difference between the value V of the
		variable at a given grid point and its value Vo at the
		observing station; ND = (V - Vo) / (V + Vo).

MAXND(0.6)
Maximum absolute value of normalized difference. Normalized differences exceeding MAXND will be replaced by MAXND, with the same sign. This applies both to normalized difference values between pairs of stations used for computing the regression coefficients, and for the values used to compute output pixels. The default value, 0.6, corresponds to a maximum ratio between values of 4 to 1.

MINFRACT(0.52)
Minimum fraction of stations above threshold for non-fill output. If MINFRACT is greater than zero, the program determines the weighted fraction of observing stations contributing to each output pixel for which the observed variable has a value greater than or equal to THRESH. If this fraction is less than MINFRACT, the output pixel's value will be the value specified by FILL.

THRESH(0.001)
Threshold for using variable value from station. If the observed value of the variable at the station is below THRESH, the station will not be used for computing the weighted value of the output pixel.

FILL(0.0)
Fill value for output pixels failing MINFRACT test.

PRINT(TERM)
Destination(s) for processing report. This summarizes the average, maximum, and minimum sum of the processing weights for each grid cell, the adjusted cutoff distance from a point to a station, and the values of slope and intercept used. This listing may be sent to one or more destinations as specifed by the value(s) entered for PRINT. Options are


    TERM:	Listing is displayed on the user's terminal.
    LP:		Listing is sent to the  printer defined by $PRINTER.
    filename:	Listing is stored in a file with the specified filename;
		the file extension ".prt" is appended to the filename
		specified by the user.

Examples:

  1. LAS> surface-temp inelev="mean_elev_1km(101, 201, 200, 400)" inobs=maxt_2001_06_18.dat out=maxt_2001_06_18

    The observed maximum temperatures for June 18, 2001 at the set set of stations included in the input .DAT file are used to create a temperature surface image named MAXT_2001_06_18 for a 200 by 400 window in the region covered by the elevation grid MEAN_ELEV_1KM. A linear regression is used to adjust the temperature values for elevation, and the other standard default values for temperature-like data are also used.

  2. LAS> surface-precip inelev=smoothed_elev inobs=smoothed_precip_6.dat out=interp_precip_6 maxnd=.333

    The temporarlly smoothed precipitation in the input .DAT file is interpolated to the region covered by the spatially smoothed elevation grid specified by image SMOOTHED_ELEV. The results are written to image INTERP_PRECIP_6. The precipitaion values are adjusted for elevation using a regression of the normalized difference between precipitation at different points against the differrence in elevation. The normalized difference is constrained by the value of MAXND to lie between -.333 and plus .333; this corresponds to a maximum elevation correction factor of 2.0 for the precipitation. Temporal smoothing of the precipitation values and spatial smoothing of the elevation values are assumed to have been done in previous processing steps -- see Usernote 1.

Description/Algorithm:

This program uses a Gaussian weighting function with a spatially varying cutoff distance defined, for station i, as

w(i) = exp(-alpha * (d / Dmax) ** 2) - exp(-alpha), d <= Dmax,

w(i) = 0, d > Dmax

where w(i) is the weight for station i, d is the distance from the station to the grid point for which the interpolation is being done, Dmax is the cutoff distance for this grid point, and alpha is the shape parameter for the Guassian.

After getting the program arguments from the TAE interface, the program reads the station observations file, and uses the latitude and longitude for each station to determine line and sample numbers of the INELEV image pixel whose center is closest to the station location.

SURFACE then defines a coarse grid covering the region defined by image INELEV such that there are an average of about 25 grid points for each observing station. At each point of this coarse grid, an iterative method is used to determine the cutoff distance Dmax which causes the number of stations within distance Dmax of the grid point to be approximately equal to the target number of stations specified by program argument STACNT.

The program then creates the output image grid by computing the weighted, interpolated value of the observed variable at the center of each pixel in the input elevation image as follows: First, the cutoff distance for the output point is computed from the coarse grid of cutoff values using bilinear interpolation. This distance is used to identify all observing stations which are within the cutoff distance of the output point.

Second, the slope of the relation between the observed variable and elevation is determined, unless a constant value has been specified by argument SLOPE or in the station observations file, or by specifying REGRESS = "NONE" (i.e., constant slope of zero). If linear interpolation is requested, the program computes the slope of the best-fit regression line between observed value and elevation for all stations within the cutoff distance of the output point, weighting each station with the weighting function defined above. If normalized difference regression is requested, the regression is performed on the normalized difference Normdiff(i,j) between the observed values V(i) and V(j) at stations i and j, which is given by

Normdiff(i,j) = [V(i) - V(j)] / [V(i) + V(j)]

For the final step in the intepoloation process, the program uses the regression slope to adjust the observed value at each station for the difference in elevation between the station and the output point, and computes the weighted average of these adjusted values. When computing the adjusted value using the normalized difference regression, the elevation adjustment for any station is reduced, if necessary, to keep it within the limits specified by program argument MAXND. The same limits are also applied to each pair of stations when computing the regression line.

Nonfatal Error Messages: None

Fatal Error Messages:

  1. [surface-fatal] Fatal error encountered

    An error condition occurred which required aborting processing. This condition is described in the preceding message(s).

  2. [surface-prtopen] Error opening temporary report file for output

    The program was unable to open an output file. Check that user has write permission to the directory.

  3. [surface-getarg] Error getting program arguments

    The program was unable to get required program argument(s) from the TAE interface. This may be a software error. Please notify LAS system administrator.

  4. [surface-oneband] Must specify single input elevation band

    The INELEV argument may only specify a single band. If the INELEV image has more than one band, a spectral window containing a single band must be specified.

  5. [surface-malloc] Memory allocation error

    The program was unable to allocate memory needed for dynamically created array. This may be a software error. Pleas notify LAS system administrator.

  6. [surface-obsopen] Error opening observations file for input

    The file specified by argument INOBS could not be opened. Check that the file path is specified correctly and that the user has read access to it.

  7. [surface-obsread] Error reading observations file

    An error occurred while reading a record from the observations file. The file may have become corrupted.

  8. [surface-slope] Bad slope specification in observations file

    The station observations file contains the keyword "slope", but it is not followed by "=" and an integer or real value.

  9. [surface-badobs] Cannot interpret obs file

    The records in the station observation file are not in the current format. See the help entry for argument INOBS.

  10. [surface-elevunit] Invalid projection units for elevation image: <units>

    The data descriptor record contains the string specified by <units> which is not among those recognized by the LAS projection conversion software. Use EDITDDR to change it to a recognized string -- as of August, 2001, the recognized strings are FEET, MET, RAD, DEG, SEC, and DMS, or any strings which start with one of these in either upper or lower case.

  11. [surface-normd] Sum of observations is zero, values XXX and YYY

    The normalized difference between two observations cannot be computed because their sum is zero. Check that argument THRESH is above zero -- a normalized difference cannot be computed when one or both values are negative.

  12. [surface-zerocut] Zero cutoff value for line XXX, pixel YYY

    The cutoff distance for the specified line and pixel was zero. This might occur if many stations all have the same latitude and longitude values specified for them.

User Notes:

  1. This program does not apply spatial smmothing to elevation or observation data nor temporal smoothing to observations, although these may improve results for precipitation. A spatially smoothed elevation grid may be produced by using program FILTER on a standard elevation dataset. Since the input file of observations data includes both station elevations and the values of the observed variables, any desired smoothing can be done by the script which creates the obervations file.

  2. Note that, unlike inverse-distance interpolation algorithms, the exponential weighting function does not insure that the output value at the location of an observing station is the same as the value reported by the station. Instead, the value will be a weighted average of the values for all stations within the cutoff distance of the output grid cell.

  3. Although this program was specifically intended for use with elevation-dependent observations, it can potentially be applied to observations which are dependent on other surface properties. In such a case, the INELEV image would specify the particular surface property at each grid point.