Python Script to creat service areas and compare their coverage

[Anonymous User]

Original User: frankyy72

Hello Members,

am stuck on my code for service area analysis. I have a .shp of hospital locations that have been geocoded, a networked dataset that I built and a county layer for one state. The problem is, I would like to get a python script that can run the service area tool and calculate three diffrent sets of service areas around the hospitals and after that, compare the coverages of the three diffrent service area rings in terms of the number of counties covered.
Up to date, I built a model in model builder to draw the service areas and I exported it to python, but I cannot get the script to do all the things I need it to do, and I have gone through many of the scripts on service areas but it doesnt seem to do what i want. Probably I am missing many things. I am new in programming and I am just stuck. I would appreciate an help, to get me going.
I have attached a copy of the coding and also a screen pic of the model builder and the output i have currently. Thanks for helping.

[Anonymous User]

Original User: JSkinn3

The problem is, I would like to get a python script that can run the service area tool and calculate three diffrent sets of service areas around the hospitals and after that, compare the coverages of the three diffrent service area rings in terms of the number of counties covered.

Hi Frank,

How do you want this returned?  For example, do you want the service area feature class to have a field that shows a count of how many counties it covers?

[FRANKYEBOAH]

Hi Frank,

How do you want this returned?  For example, do you want the service area feature class to have a field that shows a count of how many counties it covers?

Hi Jake,

Thanks for the reply. Yes, I would like to have a field show the number of counties covered by each of the service areas, for instance the 20km travel, 30km travel etc and also if I can have a table calculate percentages of the number of fields covered, I would appreciate that. Is it also possible to get a buffreing ring to compare to the service areas in the same manner?

Thanks

[Anonymous User]

Original User: JSkinn3

Here is an example on how to get the number of counties each service area intersects:

import arcpy
from arcpy import env
env.workspace = r"C:\temp\test.gdb"
env.overwriteOutput = 1

serviceAreas = "ServiceAreas"
counties = "Counties"

arcpy.AddField_management(serviceAreas, "COUNT", "LONG")

with arcpy.da.UpdateCursor(serviceAreas, ["OBJECTID", "COUNT"]) as cursor:
   for row in cursor:
      OID = row[0]
      arcpy.MakeFeatureLayer_management(serviceAreas, "SA_lyr", "OBJECTID = " + str(OID))
      arcpy.MakeFeatureLayer_management(counties, "County_lyr")
      arcpy.SelectLayerByLocation_management("County_lyr", "INTERSECT", "SA_lyr", "#", "NEW_SELECTION")
      result = arcpy.GetCount_management("County_lyr")
      count = int(result.getOutput(0))
      row[1] = count
      cursor.updateRow(row)

I'm not sure what you mean by 'calculate percentage of the number of fields covered' or 'get a buffreing ring to compare to the service areas in the same manner'.  Could you elaborate?

[FRANKYEBOAH]

Here is an example on how to get the number of counties each service area intersects:


I'm not sure what you mean by 'calculate percentage of the number of fields covered' or 'get a buffreing ring to compare to the service areas in the same manner'.  Could you elaborate?

Hi Jake,
I meant that after getting the number of counties covered by each of the rings, can we calculate percentages? for instance if all the 20km rings covers a total of 50 counties, what percentage is that out of the total number of counties?
On the buffering, I wanted to also compare 20km buffer, 3okm buffer etc for the same counties and then finally compare it with the coverage of the service areas and see the differences between the two. In the end, i would have buffers and service areas at the same distance interval and their coverages. Thanks

[Anonymous User]

Original User: JSkinn3

To find the percentage, you will just need to add a couple more lines of code:

import arcpy
from arcpy import env
env.overwriteOutput = 1
env.workspace = r"C:\temp\test.gdb"

serviceAreas = "ServiceAreas"
counties = "Counties"

arcpy.AddField_management(serviceAreas, "COUNT", "LONG")
arcpy.AddField_management(serviceAreas, "PERCENTAGE", "DOUBLE")

result = arcpy.GetCount_management(counties)
countyCount = int(result.getOutput(0))

with arcpy.da.UpdateCursor(serviceAreas, ["OBJECTID", "COUNT", "PERCENTAGE"]) as cursor:
   for row in cursor:
      OID = row[0]
      arcpy.MakeFeatureLayer_management(serviceAreas, "SA_lyr", "OBJECTID = " + str(OID))
      arcpy.MakeFeatureLayer_management(counties, "County_lyr")
      arcpy.SelectLayerByLocation_management("County_lyr", "INTERSECT", "SA_lyr", "#", "NEW_SELECTION")
      result = arcpy.GetCount_management("County_lyr")
      count = int(result.getOutput(0))
      row[1] = count
      row[2] = count / float(countyCount) * 100
      cursor.updateRow(row)

On the buffering, I wanted to also compare 20km buffer, 3okm buffer etc for the same counties

Do you want to buffer the counties or the hospitals?  One way to do this is to buffer the county/hospital and then join this to your service area feature class.  You could then create a third field and calculate the difference by dividing the area of the service area by the buffered county/hospital.

[FRANKYEBOAH]

Do you want to buffer the counties or the hospitals?  One way to do this is to buffer the county/hospital and then join this to your service area feature class.  You could then create a third field and calculate the difference by dividing the area of the service area by the buffered county/hospital.

Hi Jake,

Calculating the buffer around the hospitals would be the best and then calculate the number of counties covered by each of the rings of buffer and then calculate their percentages and finally compare all the values I get from the buffer with the values I obtained from the service area calculations, so as to understand which of them gives a meaningful representation.
Thanks!