View source for Main Page/Stuff/Luc's AML

<pre>
/*----------------------------
/* change dir
&wo D:\SAPRP_Working_Files\Raster_Calculations\Sample_Frame\costd

/*----------------------------
/* start grid
&if [show program] <> GRID &then
    grid

/* &goto skipkill1
&if [exists pointg -grid] &then
    kill pointg
&if [exists p1 -cover] &then
    kill p1
&if [exists outg -grid] &then
    kill outg
&if [exists poly1 -cover] &then
    kill poly1
&label skipkill1

/*----------------------------
/* gridpoint
&if ^ [exists frame4_n -grid] &then
    frame4_n = setnull (..\frame4 == 0, 1)
&if ^ [exists framep -point] &then
    framep = gridpoint (frame4_n)
/* add an area and count field to the frame points
&if ^ [iteminfo framep -point needlecount -exists] &then
    arc additem framep.pat framep.pat needlecount 5 5 i
&if ^ [iteminfo framep -point cdarea -exists] &then
    arc additem framep.pat framep.pat cdarea 10 10 i


/* max distance to costdistance function
&sv maxdist 100

/*----------------------------
/* count points
&describe framep
&sv numpoints = %dsc$points%

/*----------------------------
/* go!
/*&do i = 1 &to %numpoints%
&do i = 1 &to 2
    &type %i% of %numpoints%
    &type [date -full]
    /* select the point i
    clearselect framep point
    reselect framep point framep# = %i%
    /* save the selection
    writeselect pts.sel framep point
    /* reselect from framep to create p1 
    &if ^ [exists p1 -cover] &then
      arc reselect framep p1 point pts.sel

    /* delete the selection file
    &if [exists pts.sel -file] &then
      &sv d [delete pts.sel -file]
    &if [exists pts.sex -file] &then
      &sv d [delete pts.sex -file]

    /* create a grid from the one selected point  
    &if ^ [exists pointg -point] &then
      pointg = pointgrid(p1, #, #, #, 3)
    
    /* run the costdistance
    &if ^ [exists outg -grid] &then
      outg = con(costdistance(pointg, frame4_n, #, #, %maxdist%) >= 0, 1, 0)

    /* get the cell count (assumes only one value in the grid!)
    clearselect outg.vat info
    reselect outg.vat info value = 1
    &sv cellcount = [show select outg.vat info 1 item count]
    /* cell count * 9 = area
    &sv cdarea = 9 * %cellcount%

    /* create a polygon from the costdistance grid
    &if ^ [exists poly1 -cover] &then
      poly1 = gridpoly (outg, 0.001)
    
    /* reselect needle points that overlap
    clearselect needles point
    reselect needles point overlap poly1 poly 

    /* count of selection
    &sv numneedles = [extract 1 [show select needles point]]
    
    /* write selection and area back onto frameg point
    calc framep point needlecount = %numneedles%
    calc framep point cdarea = %cdarea%

    /* kill intermediate files
/* &goto skipkill
    kill pointg
    kill p1
    kill outg
    kill poly1
&label skipkill
    &type [date -full]
&end
</pre>