Wednesday, 19 June 2013

Preparing data for Kappa or Chi-square analysis

This is the first of several posts on performing Kappa and/or chi-square statistics.  These tests are commonly used to see whether remote sensing training site data matches up with that determined via image classification.  This first script prepares for subsequent analysis.

'''
PrepareData.py

Purpose:  Prepares data for Kappa and Chi-square analysis
  Author:
    Dan Patterson
    Dept of Geography and Environmental Studies
    Carleton University, Ottawa, Canada
   
Dan_Patterson@carleton.ca
  Date created:  Jan 2010
  Last modified: June 2013
  Requires:
    obs_FC       a shapefile of observed points with a classification
    buffer_dist  a buffer distance
    expected_lyr a layer file of expected values, either a shapefile or a grid.
    output_FC    an output shapefile
'''
#Import the basic modules
import os, sys, string 

import arcpy
obs_FC = sys.argv[1]
buffer_dist = float(sys.argv[2])
expected_lyr = sys.argv[3]
output_FC = sys.argv[4]

temp_FC = os.path.dirname(sys.argv[0]) + "/temp.shp"
temp2_FC = os.path.dirname(sys.argv[0]) + "/temp2.shp"
expected_type = arcpy.Describe(expected_lyr).datasetType
#
#Check to see if the expected data layer is a raster, if it is convert it to a shapefile
#Buffer the observed layer, then perform an intersection
if expected_type == "RasterDataset":
  arcpy.AddMessage("\nConverting to temporary shapefile: " + str(expected_lyr))
  arcpy.RasterToPolygon_conversion(expected_lyr, temp2_FC, "NO_SIMPLIFY", "VALUE")
  arcpy.AddMessage("\nIntersecting layers to produce: " + str(output_FC))
  arcpy.AddMessage("\nDeleting temporary files")
  if buffer_dist > 0.0:
    arcpy.Buffer_analysis(obs_FC, temp_FC, buffer_dist, "#", "#", "NONE")
    arcpy.Intersect_analysis([temp_FC, temp2_FC], output_FC, "NO_FID")
    arcpy.Delete_management(temp_FC)
  else:
    arcpy.Intersect_analysis([obs_FC, temp2_FC], output_FC, "NO_FID")
  arcpy.Delete_management(temp2_FC)
else:
  arcpy.AddMessage("\nIntersecting layers to produce: " + str(output_FC))
  arcpy.AddMessage("\nDeleting temporary files")
  if buffer_dist > 0.0:
    arcpy.Buffer_analysis(obs_FC, temp_FC, buffer_dist, "#", "#", "NONE")
    arcpy.Intersect_analysis([temp_FC, expected_lyr], output_FC,  "NO_FID")
    arcpy.Delete_management(temp_FC)
  else:
    arcpy.Intersect_analysis([obs_FC, expected_lyr], output_FC,  "NO_FID")
#
#Add an area field and calculate area.
if buffer_dist > 0.0:
  arcpy.AddField_management(output_FC, "Area", "DOUBLE", 16, 6)
  arcpy.CalculateField_management(output_FC, "Area", "!shape.area!", "PYTHON_9.3")
arcpy.AddMessage("\nProcessing complete....\n")


No comments:

Post a Comment