Main Page/Stuff/Luc's AML
From phurvitz
/*---------------------------- /* 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