"""
ArcPy Buffer Analysis Examples
============================

This script demonstrates various buffer analysis techniques using ArcPy.
Includes basic buffers, multiple distance buffers, dissolve options, and 
advanced buffer scenarios for spatial analysis.

Requirements:
- ArcGIS Pro with Spatial Analyst extension
- Valid input feature classes
- Appropriate coordinate system (projected recommended for accurate distances)
"""

import arcpy
import os
from arcpy import env

# Set workspace and environment settings
workspace = r"C:\GIS_Projects\BufferAnalysis\BufferAnalysis.gdb"
env.workspace = workspace
env.overwriteOutput = True

# Enable spatial analyst if available
try:
    arcpy.CheckOutExtension("Spatial")
    print("Spatial Analyst extension checked out successfully")
except arcpy.ExecuteError:
    print("Spatial Analyst extension not available")

def basic_buffer_example():
    """
    Example 1: Basic buffer analysis
    Creates a simple buffer around point features
    """
    print("\n=== Basic Buffer Example ===")
    
    try:
        # Input and output parameters
        input_features = "city_points"
        output_buffer = "city_buffers_5km"
        buffer_distance = "5000 Meters"  # 5km buffer
        
        # Create buffer
        arcpy.analysis.Buffer(
            in_features=input_features,
            out_feature_class=output_buffer,
            buffer_distance_or_field=buffer_distance,
            line_side="FULL",
            line_end_type="ROUND",
            dissolve_option="NONE"
        )
        
        print(f"Basic buffer created: {output_buffer}")
        
        # Get result count
        result = arcpy.GetCount_management(output_buffer)
        print(f"Number of buffer features created: {result.getOutput(0)}")
        
    except arcpy.ExecuteError:
        print(f"Error in basic buffer: {arcpy.GetMessages()}")

def multiple_distance_buffers():
    """
    Example 2: Multiple distance buffers
    Creates buffers at multiple distances (rings)
    """
    print("\n=== Multiple Distance Buffers ===")
    
    try:
        input_features = "schools"
        output_buffer = "school_service_areas"
        
        # Multiple distances: 500m, 1km, 2km
        buffer_distances = [500, 1000, 2000]
        
        buffer_list = []
        
        for i, distance in enumerate(buffer_distances):
            temp_buffer = f"temp_buffer_{distance}m"
            
            arcpy.analysis.Buffer(
                in_features=input_features,
                out_feature_class=temp_buffer,
                buffer_distance_or_field=f"{distance} Meters",
                dissolve_option="ALL"  # Dissolve overlapping buffers
            )
            
            # Add distance field
            arcpy.management.AddField(temp_buffer, "DISTANCE", "LONG")
            arcpy.management.CalculateField(temp_buffer, "DISTANCE", distance)
            
            buffer_list.append(temp_buffer)
            print(f"Created {distance}m buffer")
        
        # Merge all buffers into one feature class
        arcpy.management.Merge(buffer_list, output_buffer)
        print(f"Multiple distance buffers merged into: {output_buffer}")
        
        # Clean up temporary buffers
        for temp_buffer in buffer_list:
            arcpy.management.Delete(temp_buffer)
            
    except arcpy.ExecuteError:
        print(f"Error in multiple distance buffers: {arcpy.GetMessages()}")

def field_based_buffer():
    """
    Example 3: Field-based buffer distances
    Uses attribute values to determine buffer distance for each feature
    """
    print("\n=== Field-Based Buffer Example ===")
    
    try:
        input_features = "industrial_sites"
        output_buffer = "contamination_zones"
        distance_field = "RISK_RADIUS"  # Field containing buffer distances
        
        # Create buffer using field values
        arcpy.analysis.Buffer(
            in_features=input_features,
            out_feature_class=output_buffer,
            buffer_distance_or_field=distance_field,
            line_side="FULL",
            line_end_type="ROUND",
            dissolve_option="NONE"
        )
        
        print(f"Field-based buffer created: {output_buffer}")
        
        # Add summary statistics
        arcpy.analysis.Statistics(
            in_table=output_buffer,
            out_table="buffer_stats",
            statistics_fields=[[distance_field, "MIN"], [distance_field, "MAX"], [distance_field, "MEAN"]]
        )
        
        print("Buffer distance statistics calculated")
        
    except arcpy.ExecuteError:
        print(f"Error in field-based buffer: {arcpy.GetMessages()}")

def dissolved_buffer_analysis():
    """
    Example 4: Dissolved buffer analysis
    Creates buffers and dissolves overlapping areas
    """
    print("\n=== Dissolved Buffer Analysis ===")
    
    try:
        input_features = "hospitals"
        buffer_distance = "3000 Meters"
        
        # Buffer with dissolve
        dissolved_buffer = "hospital_service_dissolved"
        arcpy.analysis.Buffer(
            in_features=input_features,
            out_feature_class=dissolved_buffer,
            buffer_distance_or_field=buffer_distance,
            dissolve_option="ALL"
        )
        
        # Buffer without dissolve for comparison
        individual_buffer = "hospital_service_individual"
        arcpy.analysis.Buffer(
            in_features=input_features,
            out_feature_class=individual_buffer,
            buffer_distance_or_field=buffer_distance,
            dissolve_option="NONE"
        )
        
        # Compare results
        dissolved_count = int(arcpy.GetCount_management(dissolved_buffer).getOutput(0))
        individual_count = int(arcpy.GetCount_management(individual_buffer).getOutput(0))
        
        print(f"Dissolved buffers: {dissolved_count} features")
        print(f"Individual buffers: {individual_count} features")
        print(f"Reduction: {individual_count - dissolved_count} overlapping areas merged")
        
    except arcpy.ExecuteError:
        print(f"Error in dissolved buffer analysis: {arcpy.GetMessages()}")

def buffer_intersection_analysis():
    """
    Example 5: Buffer intersection analysis
    Creates buffers and analyzes intersections with other features
    """
    print("\n=== Buffer Intersection Analysis ===")
    
    try:
        # Create buffers around pollution sources
        pollution_sources = "factories"
        pollution_buffers = "pollution_buffers"
        buffer_distance = "2000 Meters"
        
        arcpy.analysis.Buffer(
            in_features=pollution_sources,
            out_feature_class=pollution_buffers,
            buffer_distance_or_field=buffer_distance,
            dissolve_option="ALL"
        )
        
        # Find residential areas within pollution buffers
        residential_areas = "residential_zones"
        affected_residential = "affected_residential"
        
        arcpy.analysis.Intersect(
            in_features=[residential_areas, pollution_buffers],
            out_feature_class=affected_residential,
            join_attributes="ALL"
        )
        
        # Calculate affected area
        arcpy.management.AddField(affected_residential, "AREA_SQKM", "DOUBLE")
        arcpy.management.CalculateGeometryAttributes(
            affected_residential, 
            [["AREA_SQKM", "AREA"]], 
            area_unit="SQUARE_KILOMETERS"
        )
        
        # Summary statistics
        arcpy.analysis.Statistics(
            in_table=affected_residential,
            out_table="pollution_impact_stats",
            statistics_fields=[["AREA_SQKM", "SUM"]]
        )
        
        print(f"Buffer intersection analysis completed")
        print(f"Affected residential areas: {affected_residential}")
        
    except arcpy.ExecuteError:
        print(f"Error in buffer intersection analysis: {arcpy.GetMessages()}")

def line_buffer_analysis():
    """
    Example 6: Line buffer analysis
    Creates buffers around linear features (roads, rivers, etc.)
    """
    print("\n=== Line Buffer Analysis ===")
    
    try:
        input_roads = "major_highways"
        
        # Different buffer distances for different road types
        # Add buffer distance field based on road classification
        arcpy.management.AddField(input_roads, "BUFFER_DIST", "DOUBLE")
        
        # Calculate buffer distances based on road type
        with arcpy.da.UpdateCursor(input_roads, ["ROAD_TYPE", "BUFFER_DIST"]) as cursor:
            for row in cursor:
                road_type = row[0]
                if road_type == "Interstate":
                    row[1] = 100  # 100m buffer for interstates
                elif road_type == "Highway":
                    row[1] = 75   # 75m buffer for highways
                elif road_type == "Arterial":
                    row[1] = 50   # 50m buffer for arterials
                else:
                    row[1] = 25   # 25m buffer for other roads
                cursor.updateRow(row)
        
        # Create buffers using the calculated distances
        road_buffers = "road_noise_zones"
        arcpy.analysis.Buffer(
            in_features=input_roads,
            out_feature_class=road_buffers,
            buffer_distance_or_field="BUFFER_DIST",
            line_side="FULL",
            line_end_type="ROUND",
            dissolve_option="LIST",
            dissolve_field=["ROAD_TYPE"]
        )
        
        print(f"Line buffer analysis completed: {road_buffers}")
        
    except arcpy.ExecuteError:
        print(f"Error in line buffer analysis: {arcpy.GetMessages()}")

def advanced_buffer_workflow():
    """
    Example 7: Advanced buffer workflow
    Combines multiple buffer techniques for complex analysis
    """
    print("\n=== Advanced Buffer Workflow ===")
    
    try:
        # Step 1: Create graduated buffers around emergency services
        emergency_services = "emergency_services"
        service_distances = [1000, 3000, 5000]  # 1km, 3km, 5km service areas
        
        merged_buffers = []
        
        for distance in service_distances:
            buffer_name = f"emergency_buffer_{distance}m"
            
            arcpy.analysis.Buffer(
                in_features=emergency_services,
                out_feature_class=buffer_name,
                buffer_distance_or_field=f"{distance} Meters",
                dissolve_option="ALL"
            )
            
            # Add service level field
            arcpy.management.AddField(buffer_name, "SERVICE_LEVEL", "TEXT")
            if distance == 1000:
                service_level = "Excellent"
            elif distance == 3000:
                service_level = "Good"
            else:
                service_level = "Basic"
                
            arcpy.management.CalculateField(
                buffer_name, 
                "SERVICE_LEVEL", 
                f"'{service_level}'"
            )
            
            merged_buffers.append(buffer_name)
        
        # Step 2: Create service area rings (donuts)
        final_service_areas = "emergency_service_areas"
        
        # Use Erase to create rings
        # Excellent service (0-1km)
        excellent_area = merged_buffers[0]  # 1km buffer
        
        # Good service (1-3km) 
        good_area = "good_service_ring"
        arcpy.analysis.Erase(
            in_features=merged_buffers[1],  # 3km buffer
            erase_features=merged_buffers[0],  # minus 1km buffer
            out_feature_class=good_area
        )
        arcpy.management.AddField(good_area, "SERVICE_LEVEL", "TEXT")
        arcpy.management.CalculateField(good_area, "SERVICE_LEVEL", "'Good'")
        
        # Basic service (3-5km)
        basic_area = "basic_service_ring"
        arcpy.analysis.Erase(
            in_features=merged_buffers[2],  # 5km buffer
            erase_features=merged_buffers[1],  # minus 3km buffer
            out_feature_class=basic_area
        )
        arcpy.management.AddField(basic_area, "SERVICE_LEVEL", "TEXT")
        arcpy.management.CalculateField(basic_area, "SERVICE_LEVEL", "'Basic'")
        
        # Step 3: Merge all service areas
        arcpy.management.Merge(
            [excellent_area, good_area, basic_area], 
            final_service_areas
        )
        
        print(f"Advanced buffer workflow completed: {final_service_areas}")
        
        # Clean up intermediate files
        for buffer in merged_buffers:
            arcpy.management.Delete(buffer)
        arcpy.management.Delete(good_area)
        arcpy.management.Delete(basic_area)
        
    except arcpy.ExecuteError:
        print(f"Error in advanced buffer workflow: {arcpy.GetMessages()}")

def buffer_quality_control():
    """
    Example 8: Buffer quality control and validation
    Checks buffer results and provides quality metrics
    """
    print("\n=== Buffer Quality Control ===")
    
    try:
        buffer_fc = "city_buffers_5km"  # From basic example
        
        # Check if buffer feature class exists
        if not arcpy.Exists(buffer_fc):
            print(f"Buffer feature class {buffer_fc} does not exist")
            return
        
        # Get spatial reference information
        desc = arcpy.Describe(buffer_fc)
        spatial_ref = desc.spatialReference
        
        print(f"Spatial Reference: {spatial_ref.name}")
        print(f"Units: {spatial_ref.linearUnitName}")
        
        if spatial_ref.type == "Geographic":
            print("WARNING: Geographic coordinate system detected.")
            print("Consider projecting to appropriate coordinate system for accurate buffers.")
        
        # Calculate actual buffer areas and compare to expected
        arcpy.management.AddField(buffer_fc, "AREA_CALC", "DOUBLE")
        arcpy.management.CalculateGeometryAttributes(
            buffer_fc,
            [["AREA_CALC", "AREA"]],
            area_unit="SQUARE_METERS"
        )
        
        # Expected area for 5km circular buffer
        expected_area = 3.14159 * (5000 ** 2)  # π * r²
        
        # Check buffer accuracy
        with arcpy.da.SearchCursor(buffer_fc, ["AREA_CALC"]) as cursor:
            areas = [row[0] for row in cursor]
        
        avg_area = sum(areas) / len(areas)
        area_difference = abs(avg_area - expected_area)
        area_percent_diff = (area_difference / expected_area) * 100
        
        print(f"Expected buffer area: {expected_area:,.0f} sq meters")
        print(f"Average actual area: {avg_area:,.0f} sq meters")
        print(f"Difference: {area_percent_diff:.2f}%")
        
        if area_percent_diff > 5:
            print("WARNING: Significant difference in buffer areas detected")
        else:
            print("Buffer areas within acceptable range")
            
    except arcpy.ExecuteError:
        print(f"Error in buffer quality control: {arcpy.GetMessages()}")

def main():
    """
    Main function to run all buffer analysis examples
    """
    print("Starting ArcPy Buffer Analysis Examples")
    print("=" * 50)
    
    # Check ArcGIS license
    if arcpy.CheckProduct("ArcInfo") == "Available":
        print("ArcGIS Advanced license available")
    elif arcpy.CheckProduct("ArcEditor") == "Available":
        print("ArcGIS Standard license available")
    elif arcpy.CheckProduct("ArcView") == "Available":
        print("ArcGIS Basic license available")
    else:
        print("No ArcGIS license available")
        return
    
    try:
        # Run buffer analysis examples
        basic_buffer_example()
        multiple_distance_buffers()
        field_based_buffer()
        dissolved_buffer_analysis()
        buffer_intersection_analysis()
        line_buffer_analysis()
        advanced_buffer_workflow()
        buffer_quality_control()
        
        print("\n" + "=" * 50)
        print("All buffer analysis examples completed successfully!")
        
    except Exception as e:
        print(f"Error running buffer analysis examples: {str(e)}")
        print(f"ArcPy messages: {arcpy.GetMessages()}")
    
    finally:
        # Check in extensions
        try:
            arcpy.CheckInExtension("Spatial")
        except:
            pass

if __name__ == "__main__":
    main()

# Additional utility functions

def create_buffer_summary_report(buffer_fc, output_table):
    """
    Creates a summary report of buffer analysis results
    """
    try:
        # Create summary statistics
        stats_fields = [
            ["Shape_Area", "SUM"],
            ["Shape_Area", "MEAN"],
            ["Shape_Area", "MIN"],
            ["Shape_Area", "MAX"],
            ["Shape_Length", "SUM"],
            ["Shape_Length", "MEAN"]
        ]
        
        arcpy.analysis.Statistics(
            in_table=buffer_fc,
            out_table=output_table,
            statistics_fields=stats_fields
        )
        
        print(f"Buffer summary report created: {output_table}")
        
    except arcpy.ExecuteError:
        print(f"Error creating buffer summary: {arcpy.GetMessages()}")

def export_buffer_to_kml(buffer_fc, output_kml):
    """
    Exports buffer results to KML for visualization
    """
    try:
        arcpy.conversion.LayerToKML(
            layer=buffer_fc,
            out_kmz_file=output_kml,
            layer_output_scale=0,
            is_composite="NO_COMPOSITE"
        )
        
        print(f"Buffer exported to KML: {output_kml}")
        
    except arcpy.ExecuteError:
        print(f"Error exporting to KML: {arcpy.GetMessages()}")

"""
Best Practices for Buffer Analysis:

1. Coordinate Systems:
   - Use projected coordinate systems for accurate distance measurements
   - Ensure all input data uses the same coordinate system
   - Consider local projections for regional analysis

2. Buffer Distance Selection:
   - Choose appropriate distances based on analysis objectives
   - Consider real-world constraints (e.g., travel time, regulations)
   - Use field-based distances for variable buffer sizes

3. Dissolve Options:
   - Use DISSOLVE_ALL for continuous coverage analysis
   - Use DISSOLVE_LIST to group by attributes
   - Use NO_DISSOLVE to maintain individual feature identity

4. Quality Control:
   - Validate buffer areas against expected values
   - Check for topology errors in complex buffer operations
   - Verify spatial reference system accuracy

5. Performance Optimization:
   - Use appropriate spatial indexes
   - Consider data generalization for large datasets
   - Process data in manageable chunks for very large analyses

6. Documentation:
   - Record buffer parameters and assumptions
   - Document coordinate system and units used
   - Maintain metadata for reproducible analysis
"""