New approach!

src/utils/wsample.jl

@@ -15,3 +15,19 @@ function wsample_single(a, w, wsum)
     end
     return @inbounds a[j]
 end
+
+"""
+weighted sampling single implementation - v2
+
+Based on https://www.aarondefazio.com/tangentially/?p=58
+"""
+function wsample_single_2(a, w, wmax)
+	n = length(w)
+	idx = rand(1:n)
+	u = wmax * rand()
+	while u > w[idx]
+		idx = rand(1:n)
+		u = wmax * rand()
+	end
+	return idx
+end