REL

Shaded relief from elevation data


EnvironmentsPYTHON :: EASI :: MODELER
Quick linksDescription :: Parameters :: Parameter descriptions :: Details :: Example :: Algorithm :: Related

Back to top

Description


Employs a user-specified light source to convert an elevation channel to a shaded relief image. Relief values are scaled from 0 to 255. The REL function is similar in function, but ANG differs in that it produces the ANGLES of incidence between a surface and a POINT source of light.
Back to top

Parameters


rel(fili, filo, dbic, dboc, dbiw, dbow, elsz, lsrc, pxsz)

Name Type Caption Length Value range
FILI * str Input file name 1 -    
FILO * str Output file name 1 -    
DBIC * List[int] Input elevation data channel 1 - 1  
DBOC * List[int] Output shaded relief channel 1 - 1  
DBIW List[int] Raster input window 0 - 4 Xoffset, Yoffset, Xsize, Ysize
DBOW List[int] Raster output window 0 - 4 Xoffset, Yoffset, Xsize, Ysize
ELSZ List[float] Elevation step size (m) 0 - 1 0.0 - 10.0
Default: 1
LSRC * List[float] Location of light source (deg) 2 - 2  
PXSZ List[float] Pixel ground-size (m) 0 - 2  

* Required parameter
Back to top

Parameter descriptions

FILI

Specifies the name of the PCIDSK image file containing data to be "shaded".

FILO

Specifies the name of the PCIDSK file to receive the shaded relief data.

The specified file must already exist.

DBIC

Specifies the channel containing the elevation data.

DBOC

Specifies the output channel to receive the shaded relief image.

If the output file is the same as the input file, the output channel must be different from the input channel.

DBIW

Specifies the raster window (Xoffset, Yoffset, Xsize, Ysize) read from the input image. If this parameter is not specified, the entire image is processed by default.

Xoffset, Yoffset define the upper-left starting pixel coordinates of the window. Xsize is the number of pixels that define the window width. Ysize is the number of lines that define the window height.

DBOW

Specifies the raster window (Xoffset, Yoffset, Xsize, Ysize) to be output. If this parameter is not specified, the entire layer is output by default.

Xoffset, Yoffset define the upper-left starting pixel coordinates of the window. Xsize is the number of pixels that define the window width. Ysize is the number of lines that define the window height.

Note: The window offset (Xoffset, Yoffset) may be different, but the window size (Xsize, Ysize) must be the same as the specified input window size (DBIW).

ELSZ

Specifies, in meters, the elevation change corresponding to a change of one gray level. If this parameter is not specified, it defaults to 1 meter. Correct calculation of the shaded relief depends on the proper specification of this value.

LSRC

Specifies, in degrees, the azimuth and elevation angles of the diffuse light source.

The first value represents the azimuth angle. This can be any non-negative number; values greater than 359 degrees are treated as MOD(360).

The second value represents the elevation angle. This value must be between 0 and 180 degrees.

PXSZ

Specifies, in meters, the X and Y dimensions of one pixel on the ground. If this parameter is not specified, REL uses the pixel size stored in the input file by CIM, APS, or GEOSET.

Correct calculation of the shaded relief depends on the proper specification of this value.

Back to top

Details

REL uses the elevation channel stored in the input file to create a shaded relief image in the output file, based upon a user-specified light source position.

Assuming that a line is drawn connecting the point source to the top-left pixel of the image, the azimuth angle is the azimuth of this line in degrees clockwise from north (top of the image); the elevation angle is the elevation of the line in degrees from the horizontal.

Correct representation of shaded relief depends on the proper specification of the elevation increment and the ground pixel size. If unspecified, the elevation step size is set to 1.

If unspecified, the ground pixel size is established by CIM, APS, or GEOSET, and is read from the input file header by REL.

The shaded gray level at a point is calculated as the cosine of the angle between the normal vector to the surface (that is, slope and aspect) and the direction of illumination, scaled to the range of 0 to 255. If a different output range is desired, the MODEL function can be used to rescale the output channel.

All surfaces not illuminated by the light source (that is, with an angle greater than 90 degrees or a negative cosine) are set to 0.

The REL and ANG functions are similar in function, but differ on the following points: ANG produces the ANGLES of incidence between a surface and a point source of light; REL produces the amount of LIGHT reflected from a surface and a point light source.

Back to top

Example

This example calculates a shaded relief image from 16-bit elevation data from channel 10 of the file "irvine.pix" and places the output in channel 8 of the same file. Each pixel in the image represents a projected ground area of 30 m x 30 m, and each increment in the gray level of the elevation image corresponds to the default elevation change of 1 meter. The sun is in the southwest at an elevation angle of 30 degrees.

from pci.rel import rel

fili	=	'irvine.pix'	# input file
filo	=	'irvine.pix'	# output file
dbic	=	[10]	# input elevation data
dboc	=	[8]	# output shaded relief image
dbiw	=	[]	# process entire image
dbow	=	[]
elsz	=	[]	# default elevation step to 1.0
lsrc	=	[225,30]	# light source location (azimuth,elevation angle)
pxsz	=	[30,30]	# pixel size, in meters

rel( fili, filo, dbic, dboc, dbiw, dbow, elsz, lsrc, pxsz )

Run NUM to print out a portion of the input image (elevations) and the output image (shaded relief values).

 NUM     Database Image Numeric Window
 
 irvine.pix;1                  [S  6PIC   512P   512L]  
 10 [16S]: C8TO16MODEL                                 

 Offset: (  459,   49)  Size: (   10,   10)

       460   461   462   463   464   465   466   467   468   469
     ______________________________________________________________
 50|   567   577   586   598   613   623   624   618   607   590
 51|   559   568   578   598   621   635   636   634   625   609
 52|   558   568   578   596   616   630   636   640   633   619
 53|   558   569   579   594   611   625   637   646   642   630
 54|   561   573   579   588   604   619   632   642   643   637
 55|   564   577   579   583   597   613   628   638   644   644
 56|   575   586   587   590   602   615   630   642   647   645
 57|   587   595   596   598   608   618   633   646   650   647
 58|   600   608   612   616   625   634   643   648   646   636
 59|   614   622   628   635   643   650   654   651   642   626

 irvine.pix;1                  [S 12PIC   512P   512L] 
 8 [ 8U]: REL    Shaded relief from elevation channel  10         

 Offset: (  459,   49)  Size: (   10,   10)
                                       
       460   461   462   463   464   465   466   467   468   469
     ____________________________________________________________
 50|   149   147   152   159   134    78    19     0     0     0
 51|   139   142   172   194   178   114    48     0     0     0
 52|   136   137   159   186   182   146   103    41     0     0
 53|   138   135   151   176   176   164   145    83     0     0
 54|   134   122   122   161   177   169   154   115    44     0
 55|   124   106    80   129   167   171   154   125    85     0
 56|   105    89    61   103   140   155   153   119    71     0
 57|    95    73    47    78   108   128   145   122    71     8
 58|    89    76    56    65    84    88    91    72    26     0
 59|   109   100    94    99   103    93    63    24     0     0
Back to top

Algorithm

We convert the azimuth and elevation angle values for the diffuse light source in LSRC to a unit vector (lx,ly,lz):

    lx = -sin(LSRC(1)) * cos(LSRC(2))
    ly =  cos(LSRC(1)) * cos(LSRC(2))
    lz =  sin(LSRC(2))

For each pixel p = (px,py,pz) we define px,py as the pixel/row location and pz as the elevation. The immediate neighbors of p are a,b,c,d where a=(px,py-1,az), b=(px-1,py,bz), etc.

        . a .
        b p c
        . d .
We calculate the shaded relief as follows:
  1. Define two vectors, V1 and V2, across the pixel in question:

       V1 = (dx,0,dzx) = (2*PXSZ(1),0,(bz-cz)*ELSZ)   (left to right)
       V2 = (0,dy,dzy) = (0,2*PXSZ(2),(az-dz)*ELSZ)   (top to bottom)
    

    where dzx and dzy are the elevation changes between pixels, and dx and dy are the projected distances between pixels. PXSZ and ELSZ are parameters specified by the user giving pixel ground size and elevation sizes. dx and dy will be twice the actual pixel ground size because the elevation differences (or gradients) are taken over 2-pixel distances.

  2. Find the normal vector to the plane formed by the two vectors V1 and V2, by taking the cross product:

         N = V1 X V2 = (-dy*dzx, -dx*dzy, dx*dy) = (nx, ny, nz)
    

The cosine of the angle between the normal vector and the light source vector is the dot product of the vectors divided by the sum of their magnitudes:

   cos(ang) = (nx*lx + ny*ly + nz*lz) /                  
          (sqrt(nx**2+ny**2+nz**2)*sqrt(lx**2+ly**2+lz**2))

Then:

   if (cos(ang) < 0) then
      shaded relief value = 0
   else
      shaded relief value = cos(ang) * 255
   endif

© PCI Geomatics Enterprises, Inc.®, 2026. All rights reserved.