-
Notifications
You must be signed in to change notification settings - Fork 0
/
CalculatePerpendicularLinesAtLineStartNodes.py
177 lines (157 loc) · 6.39 KB
/
CalculatePerpendicularLinesAtLineStartNodes.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
try:
import sys, string, os, arcpy, traceback, math
def AddAndSubtractRadians(theta):
return (theta + 1.57079632679, theta - 1.57079632679)
def Distance(x1, y1, x2, y2):
'''Cartesian distance formula'''
return float(math.pow(((math.pow((x2-x1),2)) + (math.pow((y2 - y1),2))),.5))
def CartesianToPolar(xy1, xy2):
'''Given coordinate pairs as two lists or tuples, return the polar
coordinates with theta in radians. Values are in true radians along the
unit-circle, for example, 3.14 and not -0 like a regular python
return.'''
try:
x1, y1, x2, y2 = float(xy1[0]), float(xy1[1]), float(xy2[0]), float(xy2[1])
xdistance, ydistance = x2 - x1, y2 - y1
distance = math.pow(((math.pow((x2 - x1),2)) + (math.pow((y2 - y1),2))),.5)
if xdistance == 0:
if y2 > y1:
theta = math.pi/2
else:
theta = (3*math.pi)/2
elif ydistance == 0:
if x2 > x1:
theta = 0
else:
theta = math.pi
else:
theta = math.atan(ydistance/xdistance)
if xdistance > 0 and ydistance < 0:
theta = 2*math.pi + theta
if xdistance < 0 and ydistance > 0:
theta = math.pi + theta
if xdistance < 0 and ydistance < 0:
theta = math.pi + theta
return [distance, theta]
except:
gPrint("Error in CartesianToPolar()")
def PolarToCartesian(polarcoords):
'''A tuple, or list, of polar values(distance, theta in radians) are
converted to cartesian coords'''
r = polarcoords[0]
theta = polarcoords[1]
x = r * math.cos(theta)
y = r * math.sin(theta)
return [x, y]
def gPrint(string):
#print string
arcpy.AddMessage(string)
def stringToBoolean(x):
if x == 'true':
x = True
else:
x = False
return x
gPrint("\n Runing...Create Perpendicular Lines At The Start Point of Line Features")
gPrint(" A Two-Bit Algorithms product.\n Copyright 2011 Gerry Gabrisch (gerry@gabrisch.us)\n")
infc = arcpy.GetParameterAsText(0)
distance = arcpy.GetParameterAsText(1)
fcname = arcpy.GetParameterAsText(2)
leftside = arcpy.GetParameterAsText(3)
rightside = arcpy.GetParameterAsText(4)
bothsides = arcpy.GetParameterAsText(5)
makeperptoend = arcpy.GetParameterAsText(6)
leftside = stringToBoolean(leftside)
rightside = stringToBoolean(rightside)
bothsides = stringToBoolean(bothsides)
makeperptoend = stringToBoolean(makeperptoend)
#Get the input line files geometry as a python list.
gPrint(" Cracking Features...")
desc = arcpy.Describe(infc)
shapefieldname = desc.ShapeFieldName
rows = arcpy.SearchCursor(infc)
listofpointgeometry = []
gPrint(" Finding Line End Points...")
for row in rows:
feat = row.getValue(shapefieldname)
partnum = 0
partcount = feat.partCount
thisrecordsgeometry = []
while partnum < partcount:
part = feat.getPart(partnum)
pnt = part.next()
pntcount = 0
while pnt:
thetuple = [pnt.X, pnt.Y]
thisrecordsgeometry.append(thetuple)
pnt = part.next()
pntcount += 1
partnum += 1
startnode = [thisrecordsgeometry[0][0], thisrecordsgeometry[0][1]]
if makeperptoend:
endnode = [thisrecordsgeometry[1][0], thisrecordsgeometry[1][1]]
else:
endnode = [thisrecordsgeometry[-1][0], thisrecordsgeometry[-1][1]]
listofpointgeometry.append([startnode,endnode])
#Create the feature class to store the new geometery....
featureList = []
array = arcpy.Array()
pnt = arcpy.Point()
gPrint(" Defining Perpendicular Lines...")
for pt in listofpointgeometry:
startx = pt[0][0]
starty = pt[0][1]
endx = pt[1][0]
endy = pt[1][1]
#get a theta
polarcoor = CartesianToPolar((startx,starty), (endx,endy))
#Add and subtract the 90 degrees in radians from the line...
ends = AddAndSubtractRadians(polarcoor[1])
firstend = PolarToCartesian((float(distance),float(ends[0])))
secondend = PolarToCartesian((float(distance),float(ends[1])))
firstx2 = startx + firstend[0]
firsty2 = starty + firstend[1]
secondx2 = startx + secondend[0]
secondy2 = starty + secondend[1]
if bothsides:
pnt.X, pnt.Y = firstx2 , firsty2
array.add(pnt)
pnt.X, pnt.Y = secondx2 , secondy2
array.add(pnt)
polyline = arcpy.Polyline(array)
array.removeAll()
featureList.append(polyline)
if leftside:
pnt.X, pnt.Y = startx , starty
array.add(pnt)
pnt.X, pnt.Y = firstx2 , firsty2
array.add(pnt)
polyline = arcpy.Polyline(array)
array.removeAll()
featureList.append(polyline)
if rightside:
pnt.X, pnt.Y = pnt.X, pnt.Y = startx , starty
array.add(pnt)
pnt.X, pnt.Y = secondx2 , secondy2
array.add(pnt)
polyline = arcpy.Polyline(array)
array.removeAll()
featureList.append(polyline)
gPrint(" Creating New Feature Class...")
arcpy.CopyFeatures_management(featureList, fcname)
spatialRef = arcpy.Describe(infc).spatialReference
arcpy.DefineProjection_management(fcname, spatialRef)
gPrint (" Done")
except arcpy.ExecuteError:
msgs = arcpy.GetMessages(2)
arcpy.AddError(msgs)
gPrint(msgs)
except:
tb = sys.exc_info()[2]
tbinfo = traceback.format_tb(tb)[0]
pymsg = "PYTHON ERRORS:\nTraceback info:\n" + tbinfo + "\nError Info:\n" + str(sys.exc_info()[1])
msgs = "ArcPy ERRORS:\n" + arcpy.GetMessages(2) + "\n"
arcpy.AddError(pymsg)
arcpy.AddError(msgs)
gPrint(pymsg + "\n")
gPrint (msgs)