<s> import os , sys , math , traceback , types , datetime , time <EOL> import arcpy <EOL> from arcpy import env <EOL> from arcpy import sa <EOL> from arcpy . sa import * <EOL> from time import mktime , gmtime , strftime <EOL> timezone = { <EOL> '<STR_LIT>' : - <NUM_LIT:12> , <EOL> '<STR_LIT>' : - <NUM_LIT:11> , <EOL> '<STR_LIT>' : - <NUM_LIT:10> , <EOL> '<STR_LIT>' : - <NUM_LIT:9> , <EOL> '<STR_LIT>' : - <NUM_LIT:8> , <EOL> '<STR_LIT>' : - <NUM_LIT:7> , <EOL> '<STR_LIT>' '<STR_LIT>' : - <NUM_LIT:6> , <EOL> '<STR_LIT>' : - <NUM_LIT:5> , <EOL> '<STR_LIT>' : - <NUM_LIT:4> , <EOL> '<STR_LIT>' : - <NUM_LIT> , <EOL> '<STR_LIT>' : - <NUM_LIT:3> , <EOL> '<STR_LIT>' : - <NUM_LIT:2> , <EOL> '<STR_LIT>' : - <NUM_LIT:1> , <EOL> '<STR_LIT>' : <NUM_LIT:0> , <EOL> '<STR_LIT>' '<STR_LIT>' : <NUM_LIT:1> , <EOL> '<STR_LIT>' '<STR_LIT>' : <NUM_LIT:2> , <EOL> '<STR_LIT>' : <NUM_LIT:3> , <EOL> '<STR_LIT>' : <NUM_LIT> , <EOL> '<STR_LIT>' : <NUM_LIT:4> , <EOL> '<STR_LIT>' : <NUM_LIT> , <EOL> '<STR_LIT>' : <NUM_LIT:5> , <EOL> '<STR_LIT>' : <NUM_LIT> , <EOL> '<STR_LIT>' : <NUM_LIT:6> , <EOL> '<STR_LIT>' : <NUM_LIT:7> , <EOL> '<STR_LIT>' '<STR_LIT>' : <NUM_LIT:8> , <EOL> '<STR_LIT>' : <NUM_LIT:9> , <EOL> '<STR_LIT>' : <NUM_LIT> , <EOL> '<STR_LIT>' : <NUM_LIT:10> , <EOL> '<STR_LIT>' : <NUM_LIT:11> , <EOL> '<STR_LIT>' : <NUM_LIT:12> <EOL> } <EOL> class MissingParameterError ( Exception ) : <EOL> def __init__ ( self , value ) : <EOL> self . value = value <EOL> def __str__ ( self ) : <EOL> return repr ( self . value ) <EOL> class SunPosition : <EOL> def __init__ ( self , date , point , azimuth , altitude ) : <EOL> self . date = date <EOL> self . point = point <EOL> self . azimuth = azimuth <EOL> self . altitude = altitude <EOL> def CenterPoint ( fc ) : <EOL> if debug == True : <EOL> arcpy . AddMessage ( "<STR_LIT>" + fc ) <EOL> aoiCenter = os . path . join ( env . scratchGDB , "<STR_LIT>" ) <EOL> deleteme . append ( aoiCenter ) <EOL> arcpy . FeatureToPoint_management ( fc , aoiCenter , "<STR_LIT>" ) <EOL> shapeName = arcpy . Describe ( aoiCenter ) . shapeFieldName <EOL> if debug == True : <EOL> arcpy . AddMessage ( "<STR_LIT>" + shapeName ) <EOL> rows = arcpy . SearchCursor ( aoiCenter , "<STR_LIT>" , GCS_WGS_1984 , shapeName ) <EOL> row = next ( rows ) <EOL> geom = row . getValue ( shapeName ) <EOL> point = geom . getPart ( ) <EOL> if debug == True : <EOL> arcpy . AddMessage ( "<STR_LIT>" + str ( point . X ) + "<STR_LIT:U+002CU+0020>" + str ( point . Y ) ) <EOL> return point <EOL> def date_to_julian_day ( my_date ) : <EOL> a = ( <NUM_LIT> - my_date . month ) // <NUM_LIT:12> <EOL> y = my_date . year + <NUM_LIT> - a <EOL> m = my_date . month + <NUM_LIT:12> * a - <NUM_LIT:3> <EOL> return my_date . day + ( ( <NUM_LIT> * m + <NUM_LIT:2> ) // <NUM_LIT:5> ) + <NUM_LIT> * y + y // <NUM_LIT:4> - y // <NUM_LIT:100> + y // <NUM_LIT> - <NUM_LIT> <EOL> def CalculateSunPosition ( date , observerPoint ) : <EOL> if debug == True : <EOL> arcpy . AddMessage ( "<STR_LIT>" + str ( date ) ) <EOL> twoPi = math . pi * <NUM_LIT:2> <EOL> lat = observerPoint . Y <EOL> lon = observerPoint . X <EOL> if debug == True : <EOL> arcpy . AddMessage ( "<STR_LIT>" + str ( lat ) ) <EOL> arcpy . AddMessage ( "<STR_LIT>" + str ( lon ) ) <EOL> timetuple = date . timetuple ( ) <EOL> tempHour = timetuple . tm_hour + ( timetuple . tm_min / <NUM_LIT> ) + ( timetuple . tm_sec / <NUM_LIT> ) <EOL> julianDate = date_to_julian_day ( date ) <EOL> time = julianDate - <NUM_LIT> <EOL> mnlong = <NUM_LIT> + <NUM_LIT> * time <EOL> mnlong = mnlong % <NUM_LIT> <EOL> if mnlong < <NUM_LIT:0> : <EOL> mnlong += <NUM_LIT> <EOL> mnanom = <NUM_LIT> + <NUM_LIT> * time <EOL> mnanom = mnanom % <NUM_LIT> <EOL> if mnanom < <NUM_LIT:0> : <EOL> mnanom += <NUM_LIT> <EOL> mnanom = math . radians ( mnanom ) <EOL> eclong = mnlong + <NUM_LIT> * math . sin ( mnanom ) + <NUM_LIT> * math . sin ( <NUM_LIT:2> * mnanom ) <EOL> eclong = eclong % <NUM_LIT> <EOL> if eclong < <NUM_LIT:0> : <EOL> eclong += <NUM_LIT> <EOL> oblqec = <NUM_LIT> - <NUM_LIT> * time <EOL> eclong = math . radians ( eclong ) <EOL> oblqec = math . radians ( oblqec ) <EOL> num = math . cos ( oblqec ) * math . sin ( eclong ) <EOL> den = math . cos ( eclong ) <EOL> ra = math . atan ( num / den ) <EOL> if den < <NUM_LIT:0> : <EOL> ra += math . pi <EOL> if den >= <NUM_LIT:0> and num < <NUM_LIT:0> : <EOL> ra += twoPi <EOL> dec = math . asin ( math . sin ( oblqec ) * math . sin ( eclong ) ) <EOL> gmst = <NUM_LIT> + <NUM_LIT> * time + tempHour <EOL> gmst = gmst % <NUM_LIT> <EOL> if gmst < <NUM_LIT:0> : <EOL> gmst += <NUM_LIT> <EOL> lmst = gmst + lon / <NUM_LIT:15> <EOL> lmst = lmst % <NUM_LIT> <EOL> if lmst < <NUM_LIT:0> : <EOL> lmst += <NUM_LIT> <EOL> lmst = lmst * <NUM_LIT> <EOL> lmst = math . radians ( lmst ) <EOL> ha = lmst - ra <EOL> if ha < ( - math . pi ) : <EOL> ha += twoPi <EOL> if ha > math . pi : <EOL> ha -= twoPi <EOL> latRad = math . radians ( lat ) <EOL> zenithAngle = math . acos ( math . sin ( latRad ) * math . sin ( dec ) + math . cos ( latRad ) * math . cos ( dec ) * math . cos ( ha ) ) <EOL> az = math . acos ( ( ( math . sin ( latRad ) * math . cos ( zenithAngle ) ) - math . sin ( dec ) ) / ( math . cos ( latRad ) * math . sin ( zenithAngle ) ) ) <EOL> el = math . asin ( math . sin ( dec ) * math . sin ( latRad ) + math . cos ( dec ) * math . cos ( latRad ) * math . cos ( ha ) ) <EOL> if debug == True : <EOL> arcpy . AddMessage ( "<STR_LIT>" + str ( dec ) ) <EOL> arcpy . AddMessage ( "<STR_LIT>" + str ( latRad ) ) <EOL> arcpy . AddMessage ( "<STR_LIT>" + str ( ha ) ) <EOL> cosAzPos = ( <NUM_LIT:0> <= ( math . sin ( dec ) - math . sin ( el ) * math . sin ( lat ) ) ) <EOL> sinAzNeg = ( math . sin ( az ) < <NUM_LIT:0> ) <EOL> if debug == True : <EOL> arcpy . AddMessage ( "<STR_LIT>" + str ( el ) ) <EOL> if debug == True : <EOL> arcpy . AddMessage ( "<STR_LIT>" + str ( az ) ) <EOL> el = math . degrees ( el ) <EOL> az = math . degrees ( az ) <EOL> if ha > <NUM_LIT:0> : <EOL> az = az + <NUM_LIT> <EOL> else : <EOL> az = <NUM_LIT> - az <EOL> az = az % <NUM_LIT> <EOL> if debug == True : <EOL> arcpy . AddMessage ( "<STR_LIT>" + str ( el ) ) <EOL> if debug == True : <EOL> arcpy . AddMessage ( "<STR_LIT>" + str ( az ) ) <EOL> sp = SunPosition ( date , observerPoint , az , el ) <EOL> return sp <EOL> arcpy . CheckOutExtension ( "<STR_LIT>" ) <EOL> deleteme = [ ] <EOL> debug = False <EOL> sr = arcpy . SpatialReference ( ) <EOL> sr . factoryCode = <NUM_LIT> <EOL> sr . create ( ) <EOL> GCS_WGS_1984 = sr <EOL> inputAOI = arcpy . GetParameterAsText ( <NUM_LIT:0> ) <EOL> if inputAOI == '<STR_LIT:#>' or not inputAOI : <EOL> raise MissingParameterError ( "<STR_LIT>" ) <EOL> inputElevation = arcpy . GetParameterAsText ( <NUM_LIT:1> ) <EOL> if inputElevation == '<STR_LIT:#>' or not inputElevation : <EOL> raise MissingParameterError ( "<STR_LIT>" ) <EOL> inputDatetime = arcpy . GetParameterAsText ( <NUM_LIT:2> ) <EOL> if inputDatetime == '<STR_LIT:#>' or not inputDatetime : <EOL> raise MissingParameterError ( "<STR_LIT>" ) <EOL> inputTimezone = arcpy . GetParameterAsText ( <NUM_LIT:3> ) <EOL> if inputTimezone == '<STR_LIT:#>' or not inputTimezone : <EOL> raise MissingParameterError ( "<STR_LIT>" ) <EOL> outputHillshade = arcpy . GetParameterAsText ( <NUM_LIT:4> ) <EOL> if outputHillshade == '<STR_LIT:#>' or not outputHillshade : <EOL> raise MissingParameterError ( "<STR_LIT>" ) <EOL> env . cellSize = inputElevation <EOL> env . extent = inputAOI <EOL> env . mask = inputAOI <EOL> env . snapRaster = inputElevation <EOL> env . overwriteOutput = True <EOL> env . resample = "<STR_LIT>" <EOL> env . compression = "<STR_LIT>" <EOL> env . rasterStatistics = '<STR_LIT>' <EOL> try : <EOL> timeStruct = time . strptime ( inputDatetime , "<STR_LIT>" ) <EOL> utcoffset = time . timezone <EOL> if time . daylight > <NUM_LIT:0> : <EOL> utcoffset = time . altzone <EOL> dt = datetime . datetime . fromtimestamp ( mktime ( timeStruct ) ) <EOL> if debug == True : <EOL> arcpy . AddMessage ( "<STR_LIT>" + inputTimezone + "<STR_LIT>" + str ( dt ) ) <EOL> dt = dt - datetime . timedelta ( seconds = ( timezone [ inputTimezone ] * <NUM_LIT> * <NUM_LIT> ) ) <EOL> if debug == True : <EOL> arcpy . AddMessage ( "<STR_LIT>" + str ( dt ) ) <EOL> if debug == True : <EOL> arcpy . AddMessage ( "<STR_LIT>" ) <EOL> centerPoint = CenterPoint ( inputAOI ) <EOL> if debug == True : <EOL> arcpy . AddMessage ( "<STR_LIT>" ) <EOL> sunPosition = CalculateSunPosition ( dt , centerPoint ) <EOL> if sunPosition . altitude < <NUM_LIT:0> : <EOL> if debug == True : <EOL> arcpy . AddMessage ( "<STR_LIT>" ) <EOL> outConstRaster = CreateConstantRaster ( <NUM_LIT:0> , "<STR_LIT>" ) <EOL> outConstRaster . save ( outputHillshade ) <EOL> else : <EOL> if debug == True : <EOL> arcpy . AddMessage ( "<STR_LIT>" ) <EOL> tempHillshade = arcpy . sa . Hillshade ( inputElevation , sunPosition . azimuth , sunPosition . altitude , "<STR_LIT>" ) <EOL> tempHillshade . save ( outputHillshade ) <EOL> if debug == True : <EOL> arcpy . AddMessage ( "<STR_LIT>" ) <EOL> arcpy . CalculateStatistics_management ( outputHillshade ) <EOL> except arcpy . ExecuteError : <EOL> if debug == True : <EOL> arcpy . AddMessage ( "<STR_LIT>" + str ( time . strftime ( "<STR_LIT>" , time . localtime ( ) ) ) ) <EOL> tb = sys . exc_info ( ) [ <NUM_LIT:2> ] <EOL> tbinfo = traceback . format_tb ( tb ) [ <NUM_LIT:0> ] <EOL> arcpy . AddError ( "<STR_LIT>" + tbinfo ) <EOL> msgs = arcpy . GetMessages ( ) <EOL> arcpy . AddError ( msgs ) <EOL> print ( msgs ) <EOL> except MissingParameterError as e : <EOL> if debug == True : <EOL> arcpy . AddMessage ( "<STR_LIT>" + str ( time . strftime ( "<STR_LIT>" , time . localtime ( ) ) ) ) <EOL> tb = sys . exc_info ( ) [ <NUM_LIT:2> ] <EOL> tbinfo = traceback . format_tb ( tb ) [ <NUM_LIT:0> ] <EOL> arcpy . AddError ( "<STR_LIT>" + tbinfo ) <EOL> msg = e . value + "<STR_LIT>" <EOL> arcpy . AddError ( msg ) <EOL> print ( msg ) <EOL> except : <EOL> if debug == True : <EOL> arcpy . AddMessage ( "<STR_LIT>" + str ( time . strftime ( "<STR_LIT>" , time . localtime ( ) ) ) ) <EOL> tb = sys . exc_info ( ) [ <NUM_LIT:2> ] <EOL> tbinfo = traceback . format_tb ( tb ) [ <NUM_LIT:0> ] <EOL> pymsg = "<STR_LIT>" + tbinfo + "<STR_LIT>" + str ( sys . exc_info ( ) [ <NUM_LIT:1> ] ) <EOL> msgs = "<STR_LIT>" + arcpy . GetMessages ( ) + "<STR_LIT:\n>" <EOL> arcpy . AddError ( pymsg ) <EOL> arcpy . AddError ( msgs ) <EOL> print ( ( pymsg + "<STR_LIT:\n>" ) ) <EOL> print ( msgs ) <EOL> finally : <EOL> arcpy . CheckInExtension ( "<STR_LIT>" ) <EOL> if debug == True : arcpy . AddMessage ( "<STR_LIT>" ) <EOL> for i in deleteme : <EOL> if debug == True : arcpy . AddMessage ( "<STR_LIT>" + str ( i ) ) <EOL> if arcpy . Exists ( i ) : <EOL> arcpy . Delete_management ( i ) <EOL> pass <EOL> if debug == True : arcpy . AddMessage ( "<STR_LIT>" ) </s>
