# ---------------------------------------------------------------------------
# extract_basin_prism_values.py
# Created on: 2014-07-22 18:30:25.00000 (generated by ArcGIS/ModelBuilder)
# Description: extract prism data from raster using basin shapefiles as mask
# and output a csv file for each basin
# Created by Ryan Spies (ryan.spies@amec.com)
# ---------------------------------------------------------------------------

################### User Input #####################
RFC = 'LMRFC'
basins_folder = r'P:\NWS\GIS\LMRFC\Shapefiles\basins_calib'
# if you only want to run specific basins -> list them below
# otherwise set it equal to empty list (basins_overwrite = [])
basins_overwrite = ['DOET1'] 
output_dir = 'P:\\NWS\\GIS\\' + RFC + '\\PRISM\\Model_Builder_Output\\'  # this must contain a folder for each basin (eg. FONN7)
################# End User Input ######################

# location of PRISM Raster (CONUS)
PRISM_Dataset = r'Q:\GISLibrary\PRISM\1981_2010\4km\PRISM_ppt_30yr_normal_4kmM2_annual_asc.asc'

# Import arcpy module
import arcpy
import os
import csv
arcpy.env.overwriteOutput = True

# Check out any necessary licenses
arcpy.CheckOutExtension("spatial")

# Set Geoprocessing environments
arcpy.env.scratchWorkspace = "P:\\NWS\\GIS\\Models\\10_0_tools\\Model_Output.gdb" # temporary file storage directory
arcpy.env.parallelProcessingFactor = "50"
print 'ok so far...'

#################################################################################
# find all basins in RFC or only run the specified basin list
basin_files = os.listdir(basins_folder) # list all basin shapefiles in the above specified directory
if len(basins_overwrite) != 0:
    basin_files = basins_overwrite # use the basins_overright variable to run only specified basins instead of all RFC basins
basins = []
check_dir = os.listdir(output_dir) # list all folders in output_dir
for each in basin_files:
    if each[:5] not in basins:
        basins.append(each[:5])

print 'Identified ' + str(len(basins)) + ' basins in ' + RFC + ' input directory...'        

# loop through basins
for basin in basins:
    ## Script arguments
    Basin_Boundary = basins_folder + '\\' + basin + '.shp'

    print basin 
    #Out_text = output_dir + basin + '_prism' + '.csv' 

    ## Local variables:
    Basin_Raster = 'P:\\NWS\\GIS\\Models\\10_0_tools\\Model_Output\\' + basin
    Basin_Points = 'P:\\NWS\\GIS\\Models\\10_0_tools\\Model_Output\\' + basin + '_points'
    Stats_Table = 'P:\\NWS\\GIS\\Models\\10_0_tools\\Model_Output\\prism_stats.dbf'

    ## Process: Extract by Mask
    print 'Extracting by mask...'
    arcpy.gp.ExtractByMask_sa(PRISM_Dataset, Basin_Boundary, Basin_Raster)

    ## Process: Raster to Point
    print 'Raster to point...'
    arcpy.RasterToPoint_conversion(Basin_Raster, Basin_Points, "VALUE")
    
    # Process: Summary Statistics
    print 'Calculating Summary statistics...'
    arcpy.Statistics_analysis(Basin_Points + '.shp', Stats_Table, "GRID_CODE MEAN", "")
    
    # Process: output csv file
    print 'Creating '+ basin + '_prism.csv file...'
    rows = arcpy.SearchCursor(Stats_Table)
    prism_csv = open(output_dir + basin + '_prism' + '.csv', 'wb')
    csvFile = csv.writer(prism_csv) #output csv
    fieldnames = [f.name for f in arcpy.ListFields(Stats_Table)]

    allRows = []
    for row in rows:
        rowlist = []
        for field in fieldnames:
            rowlist.append(row.getValue(field))
        allRows.append(rowlist)

    csvFile.writerow(fieldnames)
    for row in allRows:
        csvFile.writerow(row)
    row = None
    rows = None
    prism_csv.close()
    
print 'Completed grid extraction!'
