PHP Code
The PHP program used to assign treatment prescriptions is broken
into two files:
- AssignTreatments.php
- CentroidMatrix.php
The source code for AssignTreatments.php:
<?php
// Program to read a CSV file of tamarisk presence polygons, and
// assign treatments to each polygon based on attributes
//
// Ted Manahan
// NR 505 Fall 2010
require("./CentroidMatrix.php");
// The lowest average polygon cost (best accessibility).
// This is used to scale treatment costs relative to accessibility
// Future improvement: calculate this from input data
// For the Lower Gunnison River
define("BASECOST",148013.7031);
// For the John Martin SWA
//define("BASECOST",1005);
//error_reporting(E_ALL);
ini_set('display_errors', '1');
date_default_timezone_set("America/Denver");
// Some global threshold numbers to use when assigning treatments
// Aerial control is only appropriate for dense stands
$MinPctCovForAerial=75;
// Aerial control is only appropriate for large stands
// 250 acres expressed in square meters
$MinAreaForAerial=250*4046.85642;
// Mechanical only works on fairly flat ground
$MaxSlopeForMechanical=20;
// Mechanical only worthwhile for somewhat dense stands
$MinPctCovForMechanical=50;
// function to aggregate treatment groups
function AggregateTGs ($TreatmentGroups, $TamariskPolygons)
{
// Default values to use for missing data
$DefaultHeight=2;
$DefaultPct_Cov=50;
// The number of treatment groups
$TGCount=count($TreatmentGroups);
// Aggregate treatment groups
$Group=array();
for ($TG=0; $TG<$TGCount; $TG++)
{
$Group[$TG]["Area"]=0;
$Group[$TG]["SlopeMAX"]=0;
$Group[$TG]["Height"]=0;
$Group[$TG]["Pct_Cov"]=0;
$Group[$TG]["CostMEAN"]=0;
$NumPolysInTG=count($TreatmentGroups[$TG]);
for
($TGPoly=0;$TGPoly<$NumPolysInTG;$TGPoly++)
{
$Poly=$TreatmentGroups[$TG][$TGPoly];
// Get total
area of all polygons in group
$Group[$TG]["Area"]+=$TamariskPolygons[$Poly]["Area"];
// Look for
the maximum slope over all the polygons in group
if
($TamariskPolygons[$Poly]["SlopeMAX"]>
$Group[$TG]["SlopeMAX"])
$Group[$TG]["SlopeMAX"]=$TamariskPolygons[$Poly]["SlopeMAX"];
// Sum height,
calculate average later
$Height=$TamariskPolygons[$Poly]["Height"];
if
($Height<=0)$Height=$DefaultHeight;
$Group[$TG]["Height"]+=$Height;
// Sum cover
percents, calculate average later
$Pct_Cov=$TamariskPolygons[$Poly]["Pct_Cov"];
if
($Pct_Cov<=0)$Pct_Cov=$DefaultPct_Cov;
$Group[$TG]["Pct_Cov"]+=$Pct_Cov;
// Sum costs,
calculate mean later
$Group[$TG]["CostMEAN"]+=$TamariskPolygons[$Poly]["CostMEAN"];
}
// Averages over treatment group
$Group[$TG]["Height"]
=$Group[$TG]["Height"] / $NumPolysInTG;
$Group[$TG]["Pct_Cov"]
=$Group[$TG]["Pct_Cov"] / $NumPolysInTG;
$Group[$TG]["CostMEAN"]=$Group[$TG]["CostMEAN"]/ $NumPolysInTG;
}
return($Group);
}
// function to assign the treatments
function AssignTreatment($TamariskPolygons, $TreatmentGroups)
{
global $MinPctCovForAerial;
global $MinAreaForAerial;
global $MaxSlopeForMechanical;
global $MinPctCovForMechanical;
// The number of polygons and treatment groups
$NumPolygons=count($TamariskPolygons);
$TGCount=count($TreatmentGroups);
echo "AssignTreatment: Assigning treatments for
$NumPolygons polygons ".
"in $TGCount groups.\n";
// Aggregate treatment groups
//echo "AssignTreatment: Calling AggregateTGs from
AssignTreatment\n";
$GroupAttrs=AggregateTGs ($TreatmentGroups,
$TamariskPolygons);
// Assign treatment to groups
for ($GroupNo=0; $GroupNo<$TGCount; $GroupNo++) {
//echo "\tAssigning values for
treatment group $GroupNo\n";
if
(($GroupAttrs[$GroupNo]["Area"]>$MinAreaForAerial) &&
($GroupAttrs[$GroupNo]["Pct_Cov"]>$MinPctCovForAerial))
{
$GroupAttrs[$GroupNo]["Treatment"]="Aerial";
} else {
if
(($GroupAttrs[$GroupNo]["SlopeMAX"]<=$MaxSlopeForMechanical)
&&
($GroupAttrs[$GroupNo]["Pct_Cov"]>=$MinPctCovForMechanical))
{
$GroupAttrs[$GroupNo]["Treatment"]="Mechanical";
} else {
$GroupAttrs[$GroupNo]["Treatment"]="Manual";
}
}
}
// Dis-aggregate treatment groups
for ($TG=0; $TG<$TGCount; $TG++)
{
$NumPolysInTG=count($TreatmentGroups[$TG]);
for
($TGPoly=0;$TGPoly<$NumPolysInTG;$TGPoly++)
{
$Poly=$TreatmentGroups[$TG][$TGPoly];
$TamariskPolygons[$Poly]["Treatment"]=$GroupAttrs[$TG]["Treatment"];
}
}
return($TamariskPolygons);
}
function ReadTamariskPolygons($TamariskCSVFile)
{
// Read the tamarisk data file and assign it to an array
if (($handle = fopen("$TamariskCSVFile", "r")) ==
FALSE) {
die("ReadTamariskPolygons: Failed
to open file $TamariskCSVFile.");
} else {
echo "ReadTamariskPolygons:
Opened file\n";
// Get the column headers, look
for columns we want to use
$data = fgetcsv($handle, 10000,
",");
$NumColumns = count($data);
echo "ReadTamariskPolygons:
$NumColumns columns in tamarisk file\n";
for ($c=0; $c < $NumColumns;
$c++) {
$ColumnName=$data[$c];
switch
($ColumnName) {
case "VALUE":
$Value=$c;
break;
case "Height":
$Height=$c;
break;
case "Pct_Cov":
$Pct_Cov=$c;
break;
case "AREA":
$Area=$c;
break;
case "SlopeMAX":
$SlopeMAX=$c;
break;
case "SlopeMEAN":
$SlopeMean=$c;
break;
case "CostMEAN":
$CostMEAN=$c;
break;
}
}
echo "Columns:\n";
echo "\tValue=$Value,
Height=$Height, Pct_Cov=$Pct_Cov\n";
echo "\tArea=$Area,
SlopeMAX=$SlopeMAX, SlopeMean=$SlopeMean, CostMEAN=$CostMEAN\n";
// Now get the data rows and
store them in an array
$TamariskPolygons=array();
$row = 0;
while (($data = fgetcsv($handle,
1000, ",")) !== FALSE) {
for ($c=0; $c
< $NumColumns; $c++) {
//echo "Row: $row, ".$data[$c] . "<br />\n";
$AttributeValue=$data[$c];
switch ($c) {
case $Value:
$TamariskPolygons[$row]["Value"]=$AttributeValue;
break;
case $Height:
$TamariskPolygons[$row]["Height"]=$AttributeValue;
break;
case $Pct_Cov:
$TamariskPolygons[$row]["Pct_Cov"]=$AttributeValue;
break;
case $Area:
$TamariskPolygons[$row]["Area"]=$AttributeValue;
break;
case $SlopeMAX:
$TamariskPolygons[$row]["SlopeMAX"]=$AttributeValue;
break;
case $SlopeMean:
$TamariskPolygons[$row]["SlopeMean"]=$AttributeValue;
break;
case $CostMEAN:
$TamariskPolygons[$row]["CostMEAN"]=$AttributeValue;
break;
}
}
// Calculate
scaled cost, using a base of BASECOST units
$TamariskPolygons[$row]["CostMultiplier"]=
1+(($TamariskPolygons[$row]["CostMEAN"]-BASECOST)/BASECOST);
$row++;
}
fclose($handle);
}
echo "There are $row polygons.\n";
return($TamariskPolygons);
}
function ReadCosts($FileName)
// Costs are read in cost/acre. Convert to cost/square meter
// by dividing by 4046.85642
{
if (($handle = fopen("$FileName", "r")) == FALSE) {
die("ReadCosts: Failed to open
file $FileName.");
} else {
echo "ReadCosts: Opened file
$FileName\n";
$Matrix=array();
$TreatmentCount=0;
// Discard the first row - it is
the header row
$data = fgetcsv($handle, 1000,
",");
// Columns are: 0=treatment type,
1=treatment cost/acre
while (($data = fgetcsv($handle,
1000, ",")) !== FALSE) {
// Get the
treatment name and the cost
$Treatment=$data[0];
$Cost=$data[1];
$AdjustedCost=$Cost/4046.85642;
$Matrix[$Treatment]=$AdjustedCost;
echo
"\tReadCosts: Assigned treatment $Treatment cost $AdjustedCost.\n";
$TreatmentCount++;
}
fclose($handle);
$NumTreatments=$TreatmentCount;
echo "\tReadCosts: Found
$NumTreatments treatments.\n";
return($Matrix);
}
}
function PrintGroup($GroupNo, $GroupMatrix)
{
$NoPolysInGroup=count($GroupMatrix[$GroupNo]);
//echo "$GroupNo: NoPolysInGroup=$NoPolysInGroup\n";
for ($i=0;$i<$NoPolysInGroup;$i++)
{
$PolyNo=$GroupMatrix[$GroupNo][$i];
echo "$PolyNo ";
}
echo "\n";
}
function AssignGroups ($TamariskPolygons, $DistanceMatrix, $OldGroups)
// returns a new trial grouping, based on distance between polygons,
// existing grouping, and preferred treatments
{
// Aggregate treatment groups
//echo "Calling AggregateTGs from AssignGroups\n";
$GroupAttrs=AggregateTGs ($OldGroups,
$TamariskPolygons);
$NumPolygons=count($TamariskPolygons);
$OldTGCount=count($OldGroups);
$NewGroups=array();
// We are going to merge polygons into groups, and
merge groups.
// Once a group or polygon is merged, it is not
// available to be merged again. Set a flag to
indicate
// if it has been merged.
$GroupMergedFlag=array();
$PolyMergedFlag=array();
// Set a counter for each new group to indicate how
many polygons are members
// of that group. Since new groups have not yet been
assigned, this number
// starts off at zero. There can not be more new
groups than old groups.
$NewGroupPolyCount=array();
// A couter for the number of polygons in the old
groupings
$OldGroupPolyCount=array();
// Need to keep track of which groups have already
been examined
for ($i=0;$i<$OldTGCount;$i++)
{
$GroupMergedFlag[$i]=0;
$NewGroupPolyCount[$i]=0;
$OldGroupPolyCount[$i]=count($OldGroups[$i]);
}
// Also keep track of which polygons have been
examined
for ($i=0;$i<$NumPolygons;$i++)
{
$PolyMergedFlag[$i]=0;
}
// We want to consider Aerial treatment groups
first, then mechanical,
// and finally manual
$GroupOrder[0]="Aerial";
$GroupOrder[1]="Mechanical";
$GroupOrder[2]="Manual";
$NumTreatmentTypes=count($GroupOrder);
// Count of the groups in the new grouping
$NewGroupNo=0;
// need to look at each group, and figure out if it
can be joined
// with the nearest group. The processing order will
be from least expensive
// treatment type to most expensive treatment type
for
($TreatmentTypeNo=0;$TreatmentTypeNo<$NumTreatmentTypes;$TreatmentTypeNo++)
{
$TreatmentName=$GroupOrder[$TreatmentTypeNo];
for
($OldGroupNo=0;$OldGroupNo<$OldTGCount;$OldGroupNo++)
// If this group has not yet been
merged with another group,
// start a new group in the
output grouping.
// If this group has already been
merged with another group,
// do not consider it for further
grouping
{
// Only
process this group if it has the current treatment type
if
((!$GroupMergedFlag[$OldGroupNo]) &&
($TamariskPolygons[$OldGroups[$OldGroupNo][0]]["Treatment"]==$TreatmentName))
{
// Set the flag to show that this group has been
processed
$GroupMergedFlag[$OldGroupNo]=1;
// Add the polygons currently in this group to the
new grouping
// and note that these polygons have been processed.
for ($i=0; $i<$OldGroupPolyCount[$OldGroupNo];
$i++)
{
$PolyNo=$OldGroups[$OldGroupNo][$i];
$NewGroups[$NewGroupNo][$i]=$PolyNo;
$PolyMergedFlag[$PolyNo]=1;
}
$NewGroupPolyCount[$NewGroupNo]=count($NewGroups[$NewGroupNo]);
// Find the closest group to the current group
$NearestGroup=CentroidMatrix::NearestGroup
($OldGroupNo, $OldGroups,
$TamariskPolygons, $DistanceMatrix);
$SecondNearestGroup=CentroidMatrix::NthNeighborGroup
($OldGroupNo, $OldGroups,
$TamariskPolygons, $DistanceMatrix, 1);
$NearestGroupDistance=CentroidMatrix::DistanceBetweenGroups
($NearestGroup, $OldGroupNo,
$DistanceMatrix, $OldGroups);
$SecondNearestGroupDistance=CentroidMatrix::DistanceBetweenGroups
($SecondNearestGroup,
$OldGroupNo, $DistanceMatrix, $OldGroups);
/*
echo "Considering group $OldGroupNo:";
PrintGroup($OldGroupNo, $OldGroups);
echo "\tNearest
($NearestGroup:$NearestGroupDistance): ";
PrintGroup($NearestGroup, $OldGroups);
echo "\tSecond Nearest
($SecondNearestGroup:$SecondNearestGroupDistance): ";
PrintGroup($SecondNearestGroup, $OldGroups);
*/
$GroupToMerge=-1;
if ((!$GroupMergedFlag[$NearestGroup]) &&
(CanMergeGroups($OldGroupNo, $NearestGroup,
$GroupAttrs,
$NearestGroupDistance)))
{
$GroupToMerge=$NearestGroup;
echo "Merging $OldGroupNo to
nearest group $GroupToMerge\n";
} elseif ((!$GroupMergedFlag[$SecondNearestGroup])
&&
(CanMergeGroups($OldGroupNo, $SecondNearestGroup,
$GroupAttrs,
$SecondNearestGroupDistance))) {
$GroupToMerge=$SecondNearestGroup;
echo "Merging $OldGroupNo to
second nearest group $GroupToMerge\n";
}
if ($GroupToMerge>=0)
{
//echo "We can merge groups
$OldGroupNo and $GroupToMerge.\n";
// If we can merge the groups,
Add the old group to the new group
for
($g=0;$g<$OldGroupPolyCount[$GroupToMerge]; $g++)
{
$NewGroupPolyNo=$NewGroupPolyCount[$NewGroupNo];
$PolyNo=$OldGroups[$GroupToMerge][$g];
$NewGroups[$NewGroupNo][$NewGroupPolyNo]=$PolyNo;
$PolyMergedFlag[$PolyNo]=1;
$NewGroupPolyCount[$NewGroupNo]++;
}
// This group has been merged; it
is not available for
// further merges in this
iteration
$GroupMergedFlag[$GroupToMerge]=1;
} else {
// We didn't fing a group to
merge. Try to "steal" a
// Polygon from an adjacent group.
$PolyDistanceArray=CentroidMatrix::RankPolysByDistance
($OldGroupNo,
$OldGroups, $TamariskPolygons,
$DistanceMatrix, 1);
// It would be better to define a
radius to search from
// each polygon in this group,
and continue searching until
// we are "too far" away to join
polygons
// Try to steal the nearest poly
to this group, or the second
// nearest if the fist is not OK
$OKToStealPoly=0;
$PolyToSteal=$PolyDistanceArray[0]["PolyNo"];
//echo "About to try to
integrate(0) polygon $PolyToSteal into group $OldGroupNo\n";
if
((!$PolyMergedFlag[$PolyToSteal]) &&
CanMergePolyToGroup($OldGroupNo, $GroupAttrs,
$PolyToSteal, $TamariskPolygons, $DistanceMatrix))
{
$OKToStealPoly=1;
} else {
$PolyToSteal=$PolyDistanceArray[1]["PolyNo"];
//echo "About
to try to integrate (1) polygon $PolyToSteal into group $OldGroupNo\n";
if
((!$PolyMergedFlag[$PolyToSteal]) &&
CanMergePolyToGroup($OldGroupNo, $GroupAttrs,
$PolyToSteal, $TamariskPolygons,
$DistanceMatrix))
{
$OKToStealPoly=1;
}
}
if ($OKToStealPoly)
{
// Need to
remove the polygon from the old group,
$ContainingGroup=GroupContainingPoly($OldGroups, $PolyToSteal);
echo "Moving
polygon $PolyToSteal from group $ContainingGroup to group
$OldGroupNo\n";
$TempGroup=array();
$TempPolyNo=0;
$NumPolysInGroup=count($OldGroups[$ContainingGroup]);
for
($P=0;$P<$NumPolysInGroup;$P++)
{
$PolyNo=$OldGroups[$ContainingGroup][$P];
if ($PolyNo != $PolyToSteal)
{
$TempGroup[$TempPolyNo]=$PolyNo;
$TempPolyNo++;
}
}
$OldGroups[$ContainingGroup]=$TempGroup;
// It is
possible to steal the only polygon from a group
// If so,
remove the group from future consideration
if
(count($OldGroups[$ContainingGroup]==0))
$GroupMergedFlag[$ContainingGroup]=1;
// Add polygon
to the new group
$NewGroupPolyNo=$NewGroupPolyCount[$NewGroupNo];
$NewGroups[$NewGroupNo][$NewGroupPolyNo]=$PolyToSteal;
$NewGroupPolyCount[$NewGroupNo]++;
// Set the
flag to show that this polygon has been merged
$PolyMergedFlag[$PolyToSteal]=1;
}
}
$NewGroupNo++;
}
}
}
return($NewGroups);
}
function CanMergeGroups($GroupNo1, $GroupNo2, $GroupAttrs, $Distance)
// Check if the attributes for these groups are similar enough to merge
them
{
global $MinPctCovForAerial;
global $MinAreaForAerial;
global $MaxSlopeForMechanical;
global $MinPctCovForMechanical;
$CanMergeGroups=1;
if
(($GroupAttrs[$GroupNo1]["Pct_Cov"]>=$MinPctCovForAerial) &&
($GroupAttrs[$GroupNo2]["Pct_Cov"]<$MinPctCovForAerial))
$CanMergeGroups=0;
if
(($GroupAttrs[$GroupNo2]["Pct_Cov"]>=$MinPctCovForAerial) &&
($GroupAttrs[$GroupNo1]["Pct_Cov"]<$MinPctCovForAerial))
$CanMergeGroups=0;
if
(($GroupAttrs[$GroupNo1]["Pct_Cov"]>=$MinPctCovForMechanical)
&&
($GroupAttrs[$GroupNo2]["Pct_Cov"]<$MinPctCovForMechanical))
$CanMergeGroups=0;
if
(($GroupAttrs[$GroupNo2]["Pct_Cov"]>=$MinPctCovForMechanical)
&&
($GroupAttrs[$GroupNo1]["Pct_Cov"]<$MinPctCovForMechanical))
$CanMergeGroups=0;
if
(($GroupAttrs[$GroupNo1]["SlopeMAX"]>=$MaxSlopeForMechanical)
&&
($GroupAttrs[$GroupNo2]["SlopeMAX"]<$MaxSlopeForMechanical))
$CanMergeGroups=0;
if
(($GroupAttrs[$GroupNo2]["SlopeMAX"]>=$MaxSlopeForMechanical)
&&
($GroupAttrs[$GroupNo1]["SlopeMAX"]<$MaxSlopeForMechanical))
$CanMergeGroups=0;
//echo "CanMergeGroups: Groups $GroupNo1 and
$GroupNo2 merge: $CanMergeGroups.\n";
return($CanMergeGroups);
}
function CanMergePolyToGroup($GroupNo, $GroupAttrs, $PolyNo,
$TamariskPolygons, $Distance)
// Check if the attributes for this polygon is similar enough to the
group
// Currently $Distance is not used, but it may be used in the future
{
global $MinPctCovForAerial;
global $MinAreaForAerial;
global $MaxSlopeForMechanical;
global $MinPctCovForMechanical;
$CanMergeGroups=1;
if
(($GroupAttrs[$GroupNo]["Pct_Cov"]>=$MinPctCovForAerial) &&
($TamariskPolygons[$PolyNo]["Pct_Cov"]<$MinPctCovForAerial))
$CanMergeGroups=0;
if
(($TamariskPolygons[$PolyNo]["Pct_Cov"]>=$MinPctCovForAerial)
&&
($GroupAttrs[$GroupNo]["Pct_Cov"]<$MinPctCovForAerial))
$CanMergeGroups=0;
if
(($GroupAttrs[$GroupNo]["Pct_Cov"]>=$MinPctCovForMechanical)
&&
($TamariskPolygons[$PolyNo]["Pct_Cov"]<$MinPctCovForMechanical))
$CanMergeGroups=0;
if
(($TamariskPolygons[$PolyNo]["Pct_Cov"]>=$MinPctCovForMechanical)
&&
($GroupAttrs[$GroupNo]["Pct_Cov"]<$MinPctCovForMechanical))
$CanMergeGroups=0;
if
(($GroupAttrs[$GroupNo]["SlopeMAX"]>=$MaxSlopeForMechanical)
&&
($TamariskPolygons[$PolyNo]["SlopeMAX"]<$MaxSlopeForMechanical))
$CanMergeGroups=0;
if
(($TamariskPolygons[$PolyNo]["SlopeMAX"]>=$MaxSlopeForMechanical)
&&
($GroupAttrs[$GroupNo]["SlopeMAX"]<$MaxSlopeForMechanical))
$CanMergeGroups=0;
//echo "CanMergePolyToGroup: Group $GroupNo and
polygon $PolyNo merge: $CanMergeGroups.\n";
return($CanMergeGroups);
}
function PrintTreatmentGroups ($TreatmentGroups, $TamariskPolygons)
{
for ($TG=0;$TG<count($TreatmentGroups);$TG++)
{
echo "\tGroup $TG: ";
for ($TGPoly=0;
$TGPoly<count($TreatmentGroups[$TG]);$TGPoly++)
{
$Poly=$TreatmentGroups[$TG][$TGPoly];
if ($TGPoly==0)
{
$Treatment=$TamariskPolygons[$Poly]["Treatment"];
echo "$Treatment: ";
}
echo "$Poly ";
}
echo "\n";
}
}
function GroupContainingPoly($Groups, $PolyToFind)
{
//echo "GroupContainingPoly: Looking for
$PolyToFind\n";
$ContainingGroup=-1;
$GroupCount=count($Groups);
//echo "GroupContainingPoly: Number of groups to
search is $GroupCount\n";
for ($TG=0; $TG<$GroupCount; $TG++)
{
$PolyCount=count($Groups[$TG]);
//echo "GroupContainingPoly:
Group $TG has $PolyCount polys\n";
for
($Poly=0;$Poly<$PolyCount;$Poly++)
{
$PolyNo=$Groups[$TG][$Poly];
if
($PolyNo==$PolyToFind)
$ContainingGroup=$TG;
}
}
return($ContainingGroup);
}
////////////////////////////////////////////////////////////////////////////
//End of Functions
////////////////////////////////////////////////////////////////////////////
// Input files
//$TamariskCSVFile="data/NAD1983HARNColoradoSouth0503/WildlifeRefuge/Tamarisk.csv";
//$CostCSVFile="data/NAD1983HARNColoradoSouth0503/WildlifeRefuge/TreatmentCosts.csv";
//$CentroidMatrixCSVFile="data/NAD1983HARNColoradoSouth0503/WildlifeRefuge/CentroidMatrix.csv";
//$OutputFile =
"data/NAD1983HARNColoradoSouth0503/WildlifeRefuge/Prescription.csv";
$TamariskCSVFile="data/NAD1983HARNColoradoSouth0503/WesternCO/Tamarisk.csv";
$CostCSVFile="data/NAD1983HARNColoradoSouth0503/WesternCO/TreatmentCosts.csv";
$CentroidMatrixCSVFile="data/NAD1983HARNColoradoSouth0503/WesternCO/CentroidMatrix.csv";
$OutputFile =
"data/NAD1983HARNColoradoSouth0503/WesternCO/Prescription.csv";
///////////////////////////////////////
// Read in Tamarisk file
///////////////////////////////////////
$TamariskPolygons=ReadTamariskPolygons($TamariskCSVFile);
$NumPolygons=count($TamariskPolygons);
//echo "Number of polygons: $NumPolygons\n";
//print_r($TamariskPolygons);
///////////////////////////////////////
// Read in centroid matrix file
///////////////////////////////////////
$DistanceMatrix=CentroidMatrix::ReadCSVFile($CentroidMatrixCSVFile);
//print_r($DistanceMatrix);
///////////////////////////////////////
// Read in treatment cost file
///////////////////////////////////////
$TreatmentCosts=ReadCosts($CostCSVFile);
//print_r($TreatmentCosts);
///////////////////////////////////////
// Assign treatment groups
// Initially every polygon is in a separate treatment group
///////////////////////////////////////
$ProposedTreatmentGroups=array();
$TreatmentGroups=array();
for ($i=0;$i<$NumPolygons;$i++)
{
$ProposedTreatmentGroups[$i][0]=$i;
$TamariskPolygons[$i]["Group"]=$i;
}
$NewScaledTotalCost=0;
$NewTotalCost=0;
$ScaledTotalCost=1;
$TotalCost=1;
$Trials=1;
while (($NewScaledTotalCost<=$ScaledTotalCost) &&
($Trials<200)
&&
($ProposedTreatmentGroups!= $TreatmentGroups))
{
$ScaledTotalCost=$NewScaledTotalCost;
$TotalCost=$NewTotalCost;
//////////////////////////////////////
// Assign a treatment
//////////////////////////////////////
$NewTamariskPolygons=AssignTreatment($TamariskPolygons,
$ProposedTreatmentGroups);
echo "Here are the suggested Treatments:\n";
PrintTreatmentGroups ($ProposedTreatmentGroups,
$NewTamariskPolygons);
// Calculate total treatment cost
$NewTotalCost=0;
$NewScaledTotalCost=0;
for ($Poly=0;$Poly<$NumPolygons;$Poly++)
{
$Treatment=$NewTamariskPolygons[$Poly]["Treatment"];
$Area=$NewTamariskPolygons[$Poly]["Area"];
$TreatmentCost=$TreatmentCosts[$Treatment];
$CostMultiplier=$NewTamariskPolygons[$Poly]["CostMultiplier"];
$NewTotalCost+=$TreatmentCost*$Area;
$NewScaledTotalCost+=$TreatmentCost*$Area*$CostMultiplier;
}
echo "Unajusted treatment cost is $NewTotalCost.\n";
echo "Scaled treatment cost is
$NewScaledTotalCost.\n";
// Special case for first iteration
if ($ScaledTotalCost==0)
$ScaledTotalCost=$NewScaledTotalCost+1;
// If the new grouping is better than the old
grouping, try yet
// another grouping
if ($NewScaledTotalCost<=$ScaledTotalCost)
{
$TreatmentGroups=$ProposedTreatmentGroups;
$TamariskPolygons=$NewTamariskPolygons;
// Put the tamarisk polygons in
the new treatment groups
for ($TG=0;
$TG<count($TreatmentGroups);$TG++)
{
for
($Poly=0;$Poly<count($TreatmentGroups[$TG]);$Poly++)
{
$TamariskPolygons[$TreatmentGroups[$TG][$Poly]]["Group"]=$TG;
}
}
//print_r($TamariskPolygons);
$ProposedTreatmentGroups=
AssignGroups($TamariskPolygons, $DistanceMatrix, $TreatmentGroups);
//echo "Proposed treatment
groups:\n";
//print_r($ProposedTreatmentGroups);
}
$Trials++;
}
echo "The final proposed treatment groups:\n";
PrintTreatmentGroups ($TreatmentGroups, $TamariskPolygons);
echo "Unajusted treatment cost is $TotalCost.\n";
echo "Scaled treatment cost is $ScaledTotalCost.\n";
echo "Used $Trials trials to find this solution.\n";
// Write results to a csv file
$fh = fopen($OutputFile, 'w') or die("can't open file $myFile");
fwrite($fh, "VALUE,Treatment,Group\n");
for ($Polygon=0; $Polygon<$NumPolygons; $Polygon++) {
$PolyNumber=$TamariskPolygons[$Polygon]["Value"];
$Treatment =$TamariskPolygons[$Polygon]["Treatment"];
$Group =$TamariskPolygons[$Polygon]["Group"];
$StringData="$PolyNumber,$Treatment,$Group\n";
fwrite($fh, $StringData);
}
fclose($fh);
?>
The source code for CentroidMatrix.php:
<?php
// Class to manipulate a matrix of distances between polygons
//
// Ted Manahan
// NR 505 Fall 2010
//error_reporting(E_ALL);
ini_set('display_errors', '1');
// Define a maximum distance to use as an initialization value
define("MAXDIST",999999999);
class CentroidMatrix {
// function to read the .csv file and create the
matrix
function ReadCSVFile ($FileName)
{
if (($handle = fopen("$FileName",
"r")) == FALSE) {
die("Failed to
open file $FileName.");
} else {
$PrintData=0;
echo
"ReadCSVFile: Opened centroid file $FileName\n";
// How big is
the matrix
$data =
fgetcsv($handle, 1000, ",");
$num =
count($data);
echo
"ReadCSVFile: $num columns in centroid file\n";
// look for
missing columns. For some reason ArcGIS will fail to
// generate
the centroid for some polygons
$ColumnCount=0;
$HaveData=array();
for
($i=1;$i<$num;$i++)
{
if ($PrintData) echo "Line $i: ";
$CentroidNo=$data[$i];
while ($CentroidNo>$ColumnCount)
{
// Missing a centroid data point
if ($PrintData) echo "Missing
$ColumnCount ";
$HaveData[$ColumnCount]=0;
$ColumnCount++;
}
$HaveData[$ColumnCount]=1;
$ColumnCount++;
if ($PrintData) echo "$CentroidNo, \n";
}
if
($PrintData) print_r($HaveData);
$Matrix=array();
// Data starts
in the second row
$RowNumber=0;
while (($data
= fgetcsv($handle, 10000, ",")) !== FALSE) {
//print_r($data);
while (!$HaveData[$RowNumber])
{
if ($PrintData) echo "No data for
row $RowNumber: ";
for ($ColumnNumber=0;
$ColumnNumber < $ColumnCount; $ColumnNumber++)
{
if
($RowNumber==$ColumnNumber)
{
$Matrix[$RowNumber][$ColumnNumber]=0;
if ($PrintData) echo "0,";
} else {
$Matrix[$RowNumber][$ColumnNumber]=MAXDIST;
if ($PrintData) echo MAXDIST.",";
}
}
$RowNumber++;
if ($PrintData) echo "\n";
}
if ($PrintData) echo "Have data for row $RowNumber:
";
$ColumnNumber=0;
// Get the distance to each object from the current
object
for ($col=1; $col < $num; $col++) {
// Fill in missing values
while (!$HaveData[$ColumnNumber])
{
$Matrix[$RowNumber][$ColumnNumber]=MAXDIST;
if
($PrintData) echo "$ColumnNumber:999999999,";
$ColumnNumber++;
}
$AttributeValue=$data[$col];
if ($RowNumber==$ColumnNumber)
{
$Matrix[$RowNumber][$ColumnNumber]=0;
if
($PrintData) echo "$ColumnNumber:0,";
} else {
$Matrix[$RowNumber][$ColumnNumber]=$AttributeValue;
if
($PrintData) echo "$ColumnNumber:$AttributeValue,";
}
$ColumnNumber++;
}
if ($PrintData) echo "\n";
$RowNumber++;
}
fclose($handle);
//print_r($Matrix);
return($Matrix);
}
}
function NearestGroup($GroupNumber, $Groups,
$TamariskPolygons,
$DistanceMatrix)
// Return the group that is closest to GroupNumber
{
$TGCount=count($Groups);
$PolyCount=count($TamariskPolygons);
//echo "NearestGroup: Looking for
the closest group to group # $GroupNumber\n";
$NumPolysInGroup=count($Groups[$GroupNumber]);
$ClosestPolyDist=MAXDIST+1;
$CloseList=array();
// For each polygon in this
treatment group
for ($j=0;
$j<$NumPolysInGroup; $j++)
{
$CloseList[$j]=-1;
$PolyNo=$Groups[$GroupNumber][$j];
// Look for
the closest polygon
for ($k=0;
$k<$PolyCount; $k++)
{
$Distance=$DistanceMatrix[$PolyNo][$k];
if
(($TamariskPolygons[$k]["Group"]<>$GroupNumber) &&
($Distance<$ClosestPolyDist)
&&
($TamariskPolygons[$PolyNo]["Value"]<>$TamariskPolygons[$k]["Value"]))
{
$ClosestPolyDist=$Distance;
$ClosestPoly=$k;
}
}
}
$ClosestGroup=$TamariskPolygons[$ClosestPoly]["Group"];
//echo "NearestGroup: Found
$ClosestGroup\n";
return($ClosestGroup);
}
function NthNeighborGroup($GroupNumber, $Groups,
$TamariskPolygons,
$DistanceMatrix, $N=0)
// Return the group that is the Nth closest to
GroupNumber
{
$TGCount=count($Groups);
$PolyCount=count($TamariskPolygons);
//echo "NthNeighbor: Looking for
the $N neighbor to group # $GroupNumber\n";
$NumPolysInTargetGroup=count($Groups[$GroupNumber]);
$CloseList=array();
// Counter of the treatment group
we are considering
$GroupCount=0;
// For each treatment group
for ($TG=0; $TG<$TGCount;
$TG++)
{
if
($TG<>$GroupNumber)
{
//echo "Looking at TG=$TG\n";
$CloseList[$GroupCount]["GroupNo"]=-1;
$CloseList[$GroupCount]["Distance"]=MAXDIST+1;
$NumPolysInCantidateGroup=count($Groups[$TG]);
// Look at each polygon in the cantidate treatment
group
for ($k=0; $k<$NumPolysInCantidateGroup; $k++)
{
$CantidatePolyNo=$Groups[$TG][$k];
//echo
"CantidatePolyNo=$CantidatePolyNo\n";
$ClosestPolyDist=MAXDIST+1;
// and see how far it is from
each polygon in the target group
for
($P=0;$P<$NumPolysInTargetGroup; $P++)
{
$TargetPolyNo=$Groups[$GroupNumber][$P];
//echo
"TargetPolyNo=$TargetPolyNo\n";
$Distance=$DistanceMatrix[$CantidatePolyNo][$TargetPolyNo];
//echo
"Distance matrix extry $CantidatePolyNo, $TargetPolyNo is $Distance\n";
if
($Distance<$ClosestPolyDist)
{
$ClosestPolyDist=$Distance;
$ClosestPoly=$TargetPolyNo;
}
}
}
// Add the new polygon distance and group info to
the list
$j=$GroupCount;
while ($j>0 &&
$CloseList[$j-1]["Distance"]>$ClosestPolyDist)
{
$CloseList[$j]["GroupNo"]=$CloseList[$j-1]["GroupNo"];
$CloseList[$j]["Distance"]=$CloseList[$j-1]["Distance"];
$j--;
}
$CloseList[$j]["GroupNo"]=$TG;
$CloseList[$j]["Distance"]=$ClosestPolyDist;
$GroupCount++;
}
}
//print_r($CloseList);
$N_CloseGroup=$CloseList[$N]["GroupNo"];
$N_CloseDist
=$CloseList[$N]["Distance"];
//echo "NthNeighborGroup: $N
neighbor to group $GroupNumber is ".
//
"$N_CloseGroup with distance $N_CloseDist\n";
return($CloseList[$N]["GroupNo"]);
}
function PolyInGroup($PolyNo, $GroupNo, $Groups)
// Return 1 if the polygon is in this group, else 0
{
$IsInGroup=0;
$NumPolysInGroup=count($Groups[$GroupNo]);
for
($i=0;$i<$NumPolysInGroup;$i++)
{
if
($Groups[$GroupNo][$i]==$PolyNo)
$IsInGroup=1;
}
//echo "PolyInGroup: Polygon
$PolyNo in group $GroupNo: $IsInGroup\n";
return($IsInGroup);
}
function RankPolysByDistance($GroupNumber, $Groups,
$TamariskPolygons,
$DistanceMatrix, $N=0)
// Return an array of Polygons ranked by; distance
from $Groups[$GroupNumber]
{
$TGCount=count($Groups);
$PolyCount=count($TamariskPolygons);
$PolysInCurrentGroup=count($Groups[$GroupNumber]);
//echo "NthNeighborPoly: Looking
for the $N neighbor to group # $GroupNumber\n";
$CloseList=array();
// Counter for the polygon we are
considering
$PolyNumber=0;
$CloseList[$PolyNumber]["PolyNo"]=-1;
$CloseList[$PolyNumber]["Distance"]=MAXDIST+1;
// For each polygon
for ($Poly=0;
$Poly<$PolyCount; $Poly++)
{
if
(!CentroidMatrix::PolyInGroup($Poly, $GroupNumber, $Groups))
{
//echo "NthNeighborPoly: Looking at Poly=$Poly\n";
$ClosestPolyDist=MAXDIST+1;
// See how far it is from the polygon to our group
for ($P=0;$P<$PolysInCurrentGroup; $P++)
{
$TargetPolyNo=$Groups[$GroupNumber][$P];
$Distance=$DistanceMatrix[$TargetPolyNo][$Poly];
//echo "\tNthNeighborPoly:
Distance matrix extry $Poly, $TargetPolyNo is $Distance\n";
if ($Distance<$ClosestPolyDist)
{
$ClosestPolyDist=$Distance;
$ClosestPoly=$Poly;
}
}
// Add the new polygon distance and group info to
the list
$j=$PolyNumber;
while ($j>0 &&
$CloseList[$j-1]["Distance"]>$ClosestPolyDist)
{
$CloseList[$j]["PolyNo"]=$CloseList[$j-1]["PolyNo"];
$CloseList[$j]["Distance"]=$CloseList[$j-1]["Distance"];
$j--;
}
$CloseList[$j]["PolyNo"]=$ClosestPoly;
$CloseList[$j]["Distance"]=$ClosestPolyDist;
$PolyNumber++;
}
}
//print_r($CloseList);
return($CloseList);
}
function DistanceBetweenGroups ($Group1, $Group2,
$DistanceMatrix, $Groups)
{
$G1Count=count($Groups[$Group1]);
$G2Count=count($Groups[$Group2]);
$Distance=MAXDIST+1;
for ($i=0;$i<$G1Count;$i++)
{
$P1=$Groups[$Group1][$i];
for
($j=0;$j<$G2Count;$j++)
{
$P2=$Groups[$Group2][$j];
if ($DistanceMatrix[$P1][$P2]<$Distance)
$Distance=$DistanceMatrix[$P1][$P2];
}
}
return($Distance);
}
}
?>