/* ========================================================================
/* @(#)toposcale.aml   1.2   4/14/1999   10:17:56
/*-------------------------------------------------------------------------
/*            Lab for synthetic dynamic vegephenomenology    
/*-------------------------------------------------------------------------
/*   Program: TOPOSCALE.AML
/*   Purpose: This AML generates a hierarchically scaled assessment of
/*            topographic position, based on relative differences be-
/*            tween actual pixel elevation and mean local elevation at
/*            different spatial scales (search windows).
/*
/*-------------------------------------------------------------------------
/*     Usage: &r toposcale 
/*-------------------------------------------------------------------------
/*     Notes: Keep in mind that this procedure is relatively slow. When
/*            large search windows are submitted, the program gets in-
/*            creasingly slow.
/*-------------------------------------------------------------------------
/*     Input: DEM
/*    Output: unsmoothed and smoothed, hierarchically scaled topographic 
/*            exposure maps [<outfile>, and <outfile>s]
/*-------------------------------------------------------------------------
/*   History: Niklaus E Zimmermann -  9/18/1997 - Original coding
/*            Niklaus E Zimmermann -  4/14/1999 - Improved user interface
/*=========================================================================
&sys clear
&if [show program] ne GRID &then &do
    &ty
    &ty  Attention: This AML has to be started in GRID. 
    &ty             Start grid and re-run the program.
    &ty
    &goto done
  &end

&sv dem := [response 'Enter name of DEM of the site (incl. path if necessary)    ']
&if ^ [exists %dem% -grid] &then
   &return &error Grid %dem% does not exist.
&sv out := [response 'Enter name of output grid base name                        ']
&if [exists %out% -grid] &then
   &return &error Grid %out% does already exist.
&if [exists %out%s -grid] &then
   &return &error Grid %out%s does already exist.
&sv srn := [response 'Enter the search radius of the smallest window (in pixels) ']
&sv srx := [response 'Enter the search radius of the largest window (in pixels)  ']
&sv inc := [response 'Enter the increment (in pixels) to increase search radius  ']
&sv kpg := [response 'Do you want to keep the individual exposure grids (y/n)?   ']
&sv rad := %srn%
&sv cnt := 0

&do &while %rad% LE %srx%
    exp%rad% = int(%dem% - focalmean(%dem%,CIRCLE,%rad%,NODATA))
      &messages &off
    describe exp%rad%
      &messages &on
    exp%rad%nrm = int(100 * ((exp%rad% - %GRD$MEAN%) / %GRD$STDV%))
    kill exp%rad% all
    &sv rad := %rad% + %inc%
    &sv cnt := %cnt% + 1
  &end

&if %rad% gt %srx% &then &sv rad := %rad% - %inc%
&sv srd := %rad% - %inc%
&sv ctl := %cnt% - 1
&sv cth := %cnt% + 1
  &messages &off
describe exp%rad%nrm
  &messages &on
&sv xmin := %GRD$XMIN% + [calc %rad% * %GRD$DX%]
&sv ymin := %GRD$YMIN% + [calc %rad% * %GRD$DY%]
&sv xmax := %GRD$XMAX% - [calc %rad% * %GRD$DX%]
&sv ymax := %GRD$YMAX% - [calc %rad% * %GRD$DY%]
setwindow %xmin%,%ymin%,%xmax%,%ymax%
copy exp%rad%nrm tmp%cnt%
&if %kpg% EQ 'n' OR %kpg% EQ 'N' &then kill exp%rad%nrm all

&do &while %srd% GE %srn%
    if (abs(exp%srd%nrm) gt 120) tmp%ctl% = exp%srd%nrm
       else tmp%ctl% = tmp%cnt%
    endif
    &sv rad := %rad% - %inc%
    &sv srd := %srd% - %inc%
    kill tmp%cnt% all
    &if [exist tmp%cth% -grid] &then kill tmp%cth% all
    copy exp%rad%nrm tmp%cnt%
    &if %kpg% EQ 'n' OR %kpg% EQ 'N' &then kill exp%rad%nrm all
    &sv cth := %cnt%
    &sv cnt := %cnt% - 1
    &sv ctl := %ctl% - 1
 &end
 &if [exist tmp%cth% -grid] &then kill tmp%cth% all

%out%  = tmp%cnt%
%out%s = focalmean(tmp%cnt%,CIRCLE,2,NODATA)
kill tmp%cnt% all

&label done