[Closed] Find visual center of contour/polygon

it was a bug in my code… yours is right. i fixed it. thanks!

sqrt(2) is better to replace with predefined constant for better performance

I have to use this in your code

bbox = 
			(
				struct bbox 
				(
					private
						FLT_MAX = 3.40282e+038,
					
					public
						bmin = FLT_MAX * [1,1,1],
						bmax = -FLT_MAX * [1,1,1],
					
						fn size = (bmax - bmin)
				)
			)

Otherwise the FLT_MAX is undefined.

Thank you one more time.

global MapBoxOps 
(
	struct MapBoxStruct 
	(
	private
		FLT_MAX = 3.40282e+038,
		SQRT2 = sqrt 2,
	public
		
		fn getSegDistSq p a b  =
		(
			x = a.x
			y = a.y
			dx = b.x - x
			dy = b.y - y

			if (dx != 0 or dy != 0) do
			(
				t = ((p.x - x) * dx + (p.y - y) * dy) / (dx * dx + dy * dy)

				if (t > 1) then
				(
					x = b.x
					y = b.y
				)
				else if (t > 0) do 
				(
					x += dx * t
					y += dy * t
				)
			)

			dx = p.x - x
			dy = p.y - y

			dx * dx + dy * dy
		),
		
		fn pointToPolygonDist point polygon = 
		(
			inside = false
			minDistSq = FLT_MAX

			for ring in polygon do
			(
				for i = 1 to ring.count-1 do
				(
					a = ring[i]
					b = ring[i+1]
					
					if ((a.y > point.y) != (b.y > point.y) and (point.x < (b.x - a.x) * (point.y - a.y) / (b.y - a.y) + a.x)) do inside = not inside

					minDistSq = amin minDistSq (getSegDistSq point a b)
				)
			)

			(sqrt minDistSq) * (if inside then 1 else -1)
		),
				
		fn makeBBox points = 
		(
			b = box3()
			for p in points do expandToInclude b p
			b
		),

		fn makeCell center h polygon = 
		(
			d = pointToPolygonDist center polygon
			maxdist = d + h * SQRT2
			datapair center:center data:[d,h,maxdist]
		),

		fn getCentroidCell polygon =
		(
			area = 0
			c  = [0,0,0]
			ring = polygon[1]
			for i = 1 to ring.count-1 do
			(
				a = ring[i]
				b = ring[i+1]
				
				f = a.x * b.y - b.x * a.y
				c.x += (a.x + b.x) * f
				c.y += (a.y + b.y) * f
				area += f * 3
			)

			center = if area == 0 then ring[1] else c/area
			makeCell center 0 polygon
		),
		
		fn cellSorter c1 c2 = if c1.data.z < c2.data.z then -1 else if c1.data.z > c2.data.z then 1 else 0,
			
		fn polyLabel polygon precision:1 debug:false =
		(
			-- find the bounding box of the outer ring
			outer_ring = polygon[1]
			bb = makeBBox outer_ring

			size = bb.max - bb.min
			cellSize = amin size.x size.y
			if (cellSize == 0) do return bb.min
			
			h = cellSize/2.0

			-- a priority queue of cells in order of their "potential" (max distance to polygon)
			-- using cellSorter

			-- cover polygon with initial cells
			cells = #()
			
			for x = bb.min.x to bb.max.x by cellSize do
			(
				for y = bb.min.y to bb.max.y by cellSize do 
				(
					append cells (makeCell [x + h, y + h, 0] h polygon)
				)
			)
				

			-- take centroid as the first best guess
			bestCell = getCentroidCell polygon

			-- second guess: bounding box centroid
			bboxCell = makeCell (bb.min + size/2.0) 0 polygon
			if (bboxCell.data.x > bestCell.data.x) do
			(
				bestCell = bboxCell
			)

			numProbes = cells.count

			t0 = timestamp()
			h0 = heapfree
			
			qsort cells cellSorter
			while cells.count != 0 do
			(
				-- pick the most promising cell from the queue
				
				currCell = cells[cells.count]
				cells.count = cells.count - 1

				-- update the best cell if we found a better one
				if (currCell.data.x > bestCell.data.x) do
				(
					bestCell = currCell
					if (debug) do format "found best % after % probes... precision:%\n" currCell.data.x numProbes precision
				)

				-- do not drill down further if there's no chance of a better solution
				if (currCell.data.z - bestCell.data.x > precision) do
				(
					-- split the cell into four cells
					h = currCell.data.y/2.0
					append cells (makeCell [currCell.center.x - h, currCell.center.y - h, 0] h polygon)
					append cells (makeCell [currCell.center.x + h, currCell.center.y - h, 0] h polygon)
					append cells (makeCell [currCell.center.x - h, currCell.center.y + h, 0] h polygon)
					append cells (makeCell [currCell.center.x + h, currCell.center.y + h, 0] h polygon)
					numProbes += 4
					
					qsort cells cellSorter  -- !!! only here.
				)
			)

			if (debug) do
			(
				format "num probes: %\n" numProbes 
				format "best distance: %\n" bestCell.data.x
			)
			format "precision:% time:% heap:%\n" precision (timestamp() - t0) (h0 - heapfree)
			bestCell.center
		)
	)
	MapBoxOps = MapBoxStruct()
	ok
)

here is some optimization … we should do qsort only if refine cells

also I “close” ring points (set last point as the same as first)

Thank you.
Also, we can exit the while loop when the currCell.data.z - bestCell.data.x < precision, otherwise the loop will continue with the rest of the cells.

this is correct! we have to loop them all.

also I removed Cell and BBox scructs, and replaced them with built-in datapair and box3. It saves memory

Do you mean that this(fromt he original js code)

 // do not drill down further if there's no chance of a better solution
        if (cell.max - bestCell.d <= precision) continue;

means that the while loop is not terminated when this line is true?
I thought that do not drill down further if there's no chance of a better solution means that the while loop has to be stopped to avoid splitting the cells and continuing with the calculations

continue means continue the while loop… but “continue” slows down MXS code and I changed it to:

if (currCell.data.z - bestCell.data.x > precision) do < refine cells > ...

In this case why this check is performed in the while loop?
If this is true the while loop stops?
if (cell.max - bestCell.d >= precision) continue;

no…

originally was:
if d < precision continue loop without refine, if d > do refine

i changed to:

if d > precision do refine, if not skip refine and continue while

I do the same but don’t use continue because it’s slow in MXS

the comments should be changed to:
drill down if there’s a chance of a better solution