9 from __future__
import print_function
15 import matplotlib.pyplot
as plt
16 from mpl_toolkits.mplot3d
import Axes3D
19 from opals
import pyDM
20 from opals
import Segmentation, Import, Normals
23 def set_axes_radius(ax, origin, radius):
24 ax.set_xlim3d([origin[0] - radius, origin[0] + radius])
25 ax.set_ylim3d([origin[1] - radius, origin[1] + radius])
26 ax.set_zlim3d([origin[2] - radius, origin[2] + radius])
36 def set_axes_equal(ax):
44 origin = np.mean(limits, axis=1)
45 radius = 0.5 * np.max(np.abs(limits[:, 1] - limits[:, 0]))
46 set_axes_radius(ax, origin, radius)
51 def output_addinfo(info):
52 print(
" attributes: ",end=
"")
53 if info.columns() == 0:
56 for i
in range(info.columns()):
60 print(info.name(i)+
"=", end=
"")
64 print(info.get(i), end=
"")
74 def preprocessing_steps(name):
75 imp = Import.Import(inFile=name+
".laz",tileSize=100,filter=
"generic[z>265]")
78 normals = Normals.Normals(inFile=name+
".odm", normalsAlg=opals.Types.NormalsAlgorithm.simplePlane, neighbours=8, searchMode=opals.Types.SearchMode.d2, searchRadius=1)
86 filename_stem =
"strip21"
87 preprocessing_steps(filename_stem)
92 mode = opals.Types.SegmentationMethod.condClustering
95 mod = Segmentation.Segmentation()
97 filename = filename_stem+
".odm"
98 if mode == opals.Types.SegmentationMethod.planeExtraction:
102 mod.method = opals.Types.SegmentationMethod.planeExtraction
103 mod.planeExtraction.maxDist = 0.2
104 mod.planeExtraction.maxSigma = 0.15
105 mod.planeExtraction.seedCalculator =
"NormalSigma0<0.02 AND Z > 275 ? NormalSigma0 : invalid"
106 mod.filter = [
"Region[529568.8 5338773.5 529863 5338773.9 529862.6 5338642.4 529706 5338685.4 529569.1 5338753.8]"]
108 mod.sort = opals.Types.SegmentSort.descendingPSize
109 elif mode == opals.Types.SegmentationMethod.condClustering:
110 mod.inFile = filename
113 mod.method = opals.Types.SegmentationMethod.condClustering
114 mod.condClustering.criterion =
"normalX*n[0].normalX+normalY*n[0].normalY+normalZ*n[0].normalZ>0.996"
117 mod.sort = opals.Types.SegmentSort.descendingPSize
121 segManager = mod.segments
124 odm = pyDM.Datamanager.load(filename,
True,
False)
125 pi = odm.getPointIndex()
127 print(
"The segmentation process has found "+str(segManager.sizeSegments())+
" segment(s)")
131 ax = fig.add_subplot(111, projection=
'3d')
134 for seg
in segManager.getSegments():
136 print(
"Segment-Nr %d:\t %d points; COG (%.3f / %.3f / %.3f)" % (seg.getID(), seg.sizePoints(), seg.getCoG().x, seg.getCoG().y, seg.getCoG().z ) )
140 planeParams = seg.getPlaneParams()
141 assert len(planeParams) == 5
142 print(
"\tplane parameters: a=%.4f b=%.4f c=%.4f d=%.4f (sigma0=%.3f)" % tuple(planeParams) )
144 print(
"No plane params")
147 alphaShape = seg.getAlphaShape()
149 print(
"\talpha shape consists of %d parts:" % alphaShape.sizePart())
151 col = matplotlib.colors.rgb2hex(scipy.rand(3))
152 for idx, part
in enumerate(alphaShape.parts()):
153 print(
"\t\t%d part has %3d points" % (idx,part.sizePoint()), end=
"")
157 for pt
in part.points():
164 tri = ax.plot3D(vtx_x, vtx_y, vtx_z, color=col)
166 output_addinfo(part.info())
169 print(
"\touput the first %d points: " % max_pt_output)
170 ptIds = seg.getPoints()
171 for idx, id
in enumerate(ptIds):
173 print(
"\t\t%2d.point (%.3f / %.3f / %.3f)" % (idx, pt.x, pt.y, pt.z) )
174 if idx+1 >= max_pt_output:
181 plt.savefig(filename_stem+
"_alphaShapes.svg")