--- a/pkg/imports/sr.go	Wed Jul 07 18:27:40 2021 +0200
+++ b/pkg/imports/sr.go	Wed Jul 07 19:24:03 2021 +0200
@@ -642,6 +642,17 @@
 
 		feedback.Info("Building final mesh took %v.", time.Since(start))
 		feedback.Info("Store mesh.")
+	} else {
+		start = time.Now()
+		clippingPolygonBuffered.Indexify()
+		removed = make(map[int32]struct{})
+		for i, v := range xyz {
+			if clippingPolygonBuffered.IntersectionBox2D(v.Box2D()) == mesh.IntersectionOutSide {
+				removed[int32(i)] = struct{}{}
+			}
+		}
+		feedback.Info("Clipping took %v.", time.Since(start))
+		feedback.Info("Number of points to clip %d.", len(removed))
 	}
 
 	start = time.Now()
@@ -702,7 +713,7 @@
 	} else { // SurveyTypeMarking
 		if err := generateMarkingPoints(
 			ctx, tx, feedback,
-			xyz, epsg,
+			xyz, removed, epsg,
 			id,
 		); err != nil {
 			return nil, err
@@ -953,18 +964,19 @@
 	tx *sql.Tx,
 	feedback Feedback,
 	xyz mesh.MultiPointZ,
+	removed map[int32]struct{},
 	epsg uint32,
 	id int64,
 ) error {
 	log.Printf("debug: generateMarkingPoints")
 
-	min, max := mesh.MinMaxVertex(xyz.All())
+	min, max := mesh.MinMaxVertex(xyz.FilterRemoved(removed))
 
 	log.Printf("debug: min/max %.2f/%.2f\n", min.Z, max.Z)
 
 	heights := loadClassBreaks(ctx, tx, feedback, min.Z, max.Z)
 
-	classes := heights.Classify(xyz.All())
+	classes := heights.Classify(xyz.FilterRemoved(removed))
 
 	// Should not happen ... Z values over the top.
 	if n := len(classes) - 1; n > 1 && len(classes[n]) > 0 {

--- a/pkg/mesh/vertex.go	Wed Jul 07 18:27:40 2021 +0200
+++ b/pkg/mesh/vertex.go	Wed Jul 07 19:24:03 2021 +0200
@@ -163,6 +163,16 @@
 	return math.Sqrt(v.Dot(v))
 }
 
+// Box2D constructs a Box2D of this vertex.
+func (v Vertex) Box2D() Box2D {
+	return Box2D{
+		X1: v.X,
+		Y1: v.Y,
+		X2: v.X,
+		Y2: v.Y,
+	}
+}
+
 func area(a, b, c Vertex) float64 {
 	return (b.Y-a.Y)*(c.X-b.X) - (b.X-a.X)*(c.Y-b.Y)
 }