hello everyone,
i have a face in 3D and i need to find its center, i found a method but it gives only the values for a surface in 2D ?
i found this method but i need z too.
X = SUM[(Xi + Xi+1) * (Xi * Yi+1 - Xi+1 * Yi)] / 6 / area
Y = SUM[(Yi + Yi+1) * (Xi * Yi+1 - Xi+1 * Yi)] / 6 / area
how can i do that ? the face is a polygon.
thnx
From what I could find, the standard answer is to project the polygon onto the xy plane (by ignoring the z values) and use those formulas to get x and y of the center, then project the polygon onto the xz or xy plane and use the equivalent formulas to find z of the center.
By the way, the denominator should be 6*area not 6/area
but i tryed to use this formula for xy and ignoring z but the x and y calculated donât represent the x and y of the center of my polygone that i found manually.
I found these same formulas in many sources, so I think they are correct and the most likely issue (no offense intended) is that either your code or your manual calculation is in error. Could you post the code here so people can check?
One possible source of error is that the formulas require that the points are listed sequentially around the polygon. Depending on how you got them, that might not be true (SketchUp is known to reverse edges in the outer loop of a face - you have to work with the edgeuses to detect this).
Another possibility is that these formulas are for a planar polygon, even if it is placed in 3D space and not on one of the primary planes xy, xz, or yz. If your points do not all fall on a single plane, the formulas wonât work.
WARNING - THIS ONLY WORKS FOR REGULAR POLYGONS WITH EVENLY SPACED POINTS
The âcenterâ of a polygon is usually defined as the âcentroidâ which is the average of the x, y, and z co-ordinates of the bounding points. Since a polygon (by definition) is a planar surface, the resultant [x,y,z] âcenterâ will always lie on its surface. For a given array of ordered points (pti), you can iterate through them to calculate the centroid of an open or closed polygon:
def get_centroid(pti)
pt1 = pti[0]
xt = pt1.x
yt = pt1.y
zt = pt1.z
closed = 0
len = pti.length
for i in 1...len
if(pt1 != pti[i])
xt += pti[i].x
yt += pti[i].y
zt += pti[i].z
else
closed = 1
end
end
len -= closed
xt /= len
yt /= len
zt /= len
Geom::Point3d.new(xt,yt,zt)
end
Example: For a tilted square whose four points are [0,0,0], [1,0,0], [1,1,1], and [0,1,1], the average X is (0 + 1 + 1 + 0) / 4 = 0.5, the average Y is (0 + 0 + 1 + 1] / 4 = 0.5, and the average Z is (0 + 0 + 1 + 1) / 4 = 0.5 giving the centroid = [0.5,0.5,0.5].
Note: For some polygons, the centroid may lie in the plane of the polygon, yet be outside the bounding edges.
@jimhami42 , I believe you are misunderstanding the definition of âcentroidâ. It is the average of the positions of all the points in the figure. âAll the pointsâ means all the positions inside the area bounded by the figure, not simply the vertices. Mathematically it is an integral. For polygons it happens that it can be calculated based on the vertices, but not by the method you give. Read, for example, all of the Wikipedia article.
You are correct. When I woke up a bit ago, I realized my mistake and edited my answer to warn that it only works for regular polygons. @MR_ME, I hope my reply didnât add to your confusion.
@jimhami42 what do you mean by regular polygons ? POLYGONS WITH EVENLY SPACED POINTS ? as it worked for a trapeze
The following works for any SketchUp Face (I think). It has ample diagnostics included but commented out if you want to see what it is doing.
def face_centroid(face)
raise ArgumentError.new "argument must be a Face" unless face.is_a?(Sketchup::Face)
# get plane of face for later use as 2D working plane
# SketchUp returns a four-element array in which the first three elements
# are the normal unit vector of the plane, and the fourth element is
# the signed distance along the unit vector to the model origin from the
# nearest point on the plane.
plane_params = face.plane
# puts "plane_params = #{plane_params}"
plane_normal = Geom::Vector3d.new(plane_params[0..2])
plane_offset = plane_normal.clone
plane_offset.length = -plane_params[3]
# puts "plane_normal = #{plane_normal}\nplane_offset = #{plane_offset}"
# get ordered array of vertices of face
outer_loop = face.outer_loop
edgeuses = outer_loop.edgeuses
points = []
edgeuses.each do |eu|
edge = eu.edge
if(eu.reversed?)
points << edge.end.position
else
points << edge.start.position
end
end
# formulas assume the polygon repeats the start point
points << points[0]
# puts "points:\n #{points.join(",\n")}"
# calculate parameters of transform to face plane:
# origin and two orthogonal vectors to use as basis
plane_origin = ORIGIN + plane_offset
plane_x = (Y_AXIS * plane_normal).normalize
plane_y = (plane_normal * plane_x).normalize
# puts "plane_origin = #{plane_origin}\nplane_x = #{plane_x}\nplane_y = #{plane_y}"
# show the 2D coordinate axes
# ents=Sketchup.active_model.entities
# draw_x = plane_x.clone
# draw_x.length = 20
# draw_y = plane_y.clone
# draw_y.length = 20
# draw_z = plane_normal.clone
# draw_z.length = 20
# Sketchup.active_model.start_operation("axes")
# ents.add_edges(plane_origin, plane_origin + draw_x)
# ents.add_edges(plane_origin, plane_origin + draw_y)
# ents.add_edges(plane_origin, plane_origin + draw_z)
# Sketchup.active_model.commit_operation
# transform 3D vertices to 2D coordinates on face plane
plane_points = []
points.each do |point|
plane_point = [0,0]
(0..2).each do |i|
plane_point[0] += (point[i]-plane_origin[i])*plane_x[i]
plane_point[1] += (point[i]-plane_origin[i])*plane_y[i]
end
plane_points << plane_point
end
# puts "plane_points:\n #{plane_points}"
# apply 2D formulas to get x,y center in face plane
area = 0
x_sum = 0
y_sum = 0
(0...plane_points.length-1).each do |i|
# puts "i = #{i}"
area += plane_points[i][0] * plane_points[i+1][1] - plane_points[i+1][0] * plane_points[i][1]
x_sum += (plane_points[i][0] + plane_points[i+1][0])*(plane_points[i][0]*plane_points[i+1][1] - plane_points[i+1][0]*plane_points[i][1])
y_sum += (plane_points[i][1] + plane_points[i+1][1])*(plane_points[i][0]*plane_points[i+1][1] - plane_points[i+1][0]*plane_points[i][1])
end
area = area/2
center_x = x_sum/(6*area)
center_y = y_sum/(6*area)
# puts "area = #{area}\ncenter_x = #{center_x}\ncenter_y = #{center_y}"
# reverse transform 2D vertices back to 3D
x_offset = plane_x
x_offset.length=center_x
y_offset = plane_y
y_offset.length=center_y
# puts "x_offset = #{x_offset}\ny_offset = #{y_offset}"
center = plane_origin + x_offset + y_offset
# puts "center = #{center}"
end
I came up with a 3D solution that works for convex polygons ⌠the vertices must be ordered around the perimeter.
def get_centroid(pts)
total_area = 0
total_centroids = Geom::Vector3d.new(0,0,0)
third = Geom::Transformation.scaling(1.0 / 3.0)
npts = pts.length
for i in 0...(npts-2)
vec1 = Geom::Vector3d.new(pts[i+1].x - pts[0].x, pts[i+1].y - pts[0].y, pts[i+1].z - pts[0].z)
vec2 = Geom::Vector3d.new(pts[i+2].x - pts[0].x, pts[i+2].y - pts[0].y, pts[i+2].z - pts[0].z)
area = (vec1.cross vec2).length / 2.0
total_area += area
centroid = (vec1 + vec2).transform(third)
t = Geom::Transformation.scaling(area)
total_centroids += centroid.transform(t)
end
c = Geom::Transformation.scaling(1.0 / total_area)
total_centroids.transform!(c) + Geom::Vector3d.new(pts[0].x,pts[0].y,pts[0].z)
end
I havenât exhaustively tested this, but I believe it to be correct.
thnx for you answer iâll test it and get back to you
I believe it only works if the polygon is in the octant where all the coordinates have positive values. Some tests and sign corrections are needed to make it work for the general case.
A neat approach based on 3D geometry of triangles (Iâve always been a fan of vector algebra vs analytic geometry)! The convexity is so that all the cross products point the same way - otherwise the areas donât add up, some of them have to be subtracted. That could be tested by comparing the direction of each cross product with a reference one and reversing the sign of the area as required.
@MR_ME just a detail note: both of our methods assume that there are no holes (inner loops) in the face.
I ran into similar issues with my first (unpublished) attempt. I found it easier to just do the 2D transformation than to sort out all the combinations of octants.
Actually, I oriented everything about the first point to simplify things, so the value I return is relative to that ⌠adding in the first point offset seems to take care of this:
total_centroids.transform!(c) + Geom::Vector3d.new(pts[0].x,pts[0].y,pts[0].z)
I use:
def find_centroid(face)
vertices=face.vertices
vertices[vertices.length]=vertices[0] ### =loop
centroid=Geom::Point3d.new()
a_sum=0.0
x_sum=0.0
y_sum=0.0
for i in (0...vertices.length-1)
temp=(vertices[i].position.x*vertices[i+1].position.y)-(vertices[i+1].position.x*vertices[i].position.y)
a_sum+=temp
x_sum+=(vertices[i+1].position.x+vertices[i].position.x)*temp
y_sum+=(vertices[i+1].position.y+vertices[i].position.y)*temp
end###for
###area=a_sum/2
centroid.x=x_sum/(3*a_sum)
centroid.y=y_sum/(3*a_sum)
centroid.z=face.bounds.min.z
return centroid
end
Works on âflatâ faces.
Of course depending on what you want:
center = face.bounds.center
works for âregularâ facesâŚ
I added your suggestion to get this (finally):
def get_centroid(pts)
total_area = 0
total_centroids = Geom::Vector3d.new(0,0,0)
third = Geom::Transformation.scaling(1.0 / 3.0)
npts = pts.length
vec1 = Geom::Vector3d.new(pts[1].x - pts[0].x, pts[1].y - pts[0].y, pts[1].z - pts[0].z)
vec2 = Geom::Vector3d.new(pts[2].x - pts[0].x, pts[2].y - pts[0].y, pts[2].z - pts[0].z)
ref_sense = vec1.cross vec2
for i in 0...(npts-2)
vec1 = Geom::Vector3d.new(pts[i+1].x - pts[0].x, pts[i+1].y - pts[0].y, pts[i+1].z - pts[0].z)
vec2 = Geom::Vector3d.new(pts[i+2].x - pts[0].x, pts[i+2].y - pts[0].y, pts[i+2].z - pts[0].z)
vec = vec1.cross vec2
area = vec.length / 2.0
if(ref_sense.dot(vec) < 0)
area *= -1.0
end
total_area += area
centroid = (vec1 + vec2).transform(third)
t = Geom::Transformation.scaling(area)
total_centroids += centroid.transform(t)
end
c = Geom::Transformation.scaling(1.0 / total_area)
total_centroids.transform!(c) + Geom::Vector3d.new(pts[0].x,pts[0].y,pts[0].z)
end
This should work for convex and concave polygons.
A few small comments:
3D Vector math is too abstract for me. I prefer to visually see what Iâm doing during testing. So when I need to do 2D math on a face I follow this procedure:
A: create an ordered array of cpoints: one at each vertex of the face.
B: given: originpt = Geom::Point3d.new(0, 0, 0)
C: create transform: transtoxy = Geom::Transformation.new(originpt, face.normal).inverse
D: transform the cpoints (they will all now be at the same z)
E: do whatever 2D calculations are needed on the cpoint x and y positions.
F: delete the cpoints
Also remember that any true area calculation must account for faces that are partially or completely bounded by arcs. The difference is crucial in some applications.
My simple planar face solution:
EDIT: Only regular polygon face.
def avg(*args)
sum = 0
args.flatten!
args.each {|n| sum = sum + n }
sum / args.size.to_f
end
def face_ctr(face)
pts = face.outer_loop.vertices.map {|v| v.position }
Geom::Point3d::new(
avg( pts.map{|pt| pt.x }.minmax ),
avg( pts.map{|pt| pt.y }.minmax ),
avg( pts.map{|pt| pt.z }.minmax )
).project_to_plane(
face.plane # make sure point is on the face's plane
)
end
Well, there are some flies in the ointment in some of the proposed methods. For the face in the attached image, the methods return the results below the image.
> get_centroid(sel[0].outer_loop.vertices.map {|v| v.position })
(15.362717676505655, 11.646531996020338, 8.25834202291012)
> face_centroid(sel[0])
(15.362718", 11.646532", 8.258342")
> find_centroid(sel[0])
(15.362718", 11.646532", 5.084123")
> face_ctr(sel[0])
(15.156338", 13.150177", 8.190887")