Notifications
Clear all

[Closed] How to find if the object is cylinder?

hope I can contribute in plain English:

(This assumes the geo is orthogonally oriented, if it’s slanted at an angle, there’ll be problems, also there’s only verts at either ends, no middle verts)

  1. get CoM
  2. get Dist from each point to CoM
  3. if Dist are all the same it means it’s either a sphere or a cylinder (a cube is technically a 4 sided cylinder but forget that as the OP says has more than 10 sides).
  4. get BBox, here is where I need the assumption of alignment. If it’s an some angle, it’s screwed, so assume the cyl is perfectly aligned to any of the X/Y/Z axis
  5. if one dim of the BBox is longer than the other 2…bingo, it’s a cyl
  6. what if the BBox of the cyl is a perfect cube ? hmm…still thinking…

possible flaw, what if someone wants to sabotage my test by loading a lot more verts at one end of the cyl than the other ? But I suppose I can get around that by getting BBox then get its center, instead of ‘average of all verts’.

aww…shite…the centre point on the caps…that will throw me off…sorry…back to the drawing board !!

[quote=] For scaled in one axis cylidner it returns TRUE, but it should be FALSE.
Edit: the reason to return TRUE all the time is that all verts are used
to calcualte the total radius. The cylinder is still a “symetrical”
object. This is whay I use only 3 verts, fare enough, to check the
radiuses.Only open verts should be used. Valid circular cap should have all radiuses of equal length and this is another thing we need to check when calculating average radius.
If any cap radius is longer/shorter then it’s either s not a cap or not a valid cap.
If two valid caps with equal average radius are found then it’s a cylinder for sure.

Is there always two caps (one in each side)?
Can you show please your ‘FindCenterOfCircle’ function?
I think I’ve got a solution based on yours (smooth groups).

1 Reply
(@miauu)
Joined: 10 months ago

Posts: 0

The caps should be always 2. If there is only one cap the code will works the same way. My version of the script uses shortest edges, so there is no matter if those shoretest edges are part of a cap or they are part of the other side of the cylinder when there is no cap.
The “FindCenterOfCircle” function that I use is writen by Swordslayer:


fn circumcenter p1 p2 p3 =
(
   fn barycentricToWorld p1 p2 p3 u v w = (u*p1 + v*p2 + w*p3) / (u + v + w)

   local a = p3 - p2
   local b = p1 - p3
   local c = p2 - p1

   local u = (dot a a) * (dot c b)
   local v = (dot b b) * (dot c a)
   local w = (dot c c) * (dot b a)

   barycentricToWorld p1 p2 p3 u v w
)

Serejah already posted a solution with smoothing groups that splits the cylinder to 3 elements(2 caps and the rest).

This is how my code works(using Serejah’s method to find the caps or using mine). I use 3 verts(from the cap) to find the center of the circle that they forms. Then I use another vert(as far as possible) and if the distance between those vert and the found center is the same as the distance between the first 3 verts and the found center the object is marked as a cylinder, otherwise – not a cylinder. Of course I can use 2,3…n verts to copmare the distance of the initial 3 verts used to find the center of the radius to increase the chance to detect if the objecet is cylinder or not.

If there are two caps, we should be able to split into three elements. If only one, into two ones.
A cap will always be the element with less number of verts. Perpendicularity is insured by smooth of 89.9 degrees.
Just find the center of a cap with 3 points, its radius, and check for the distance of the other cap points against radius (if there aren’t extrange possibilities, just a check with a fourth point should be enougth).

Try this. Don’t really know if it has a performance improvement.
It test all cap verts for radius. If you want to decrease this number, I will take at least 8 (but think about a prism with a star-shape base. Only 1 check is not enougth).
Using PolyTools algorithm for getting elements may improve it a lot (tell me if you don’t find it).


   fn FindCenterOfCircle p1 p2 p3 =
   (
      local a = p1-p3
      local b = p2-p3
      
      lenA = length a
      lenB = length b
      crossAB = cross a b
      lenAB = length crossAB
      
      --radius = (lenA * lenB * length (a-b)) / (2 * lenAB)
      center = (cross (lenA^2 * b - lenB^2 * a) crossAB) / (2 * lenAB^2) + p3
   )
   
   fn isCylinder theNode =
   (
      if not isKindOf theNode GeometryClass do (return false)
      o = snapshotAsMesh theNode
      meshop.autoSmooth o #{1..o.numfaces} 89.9
      smoothGroups = #()
      for i=1 to o.numfaces do
      (
         sg = getFaceSmoothGroup o i
         if smoothGroups[sg] == undefined then smoothGroups[sg] = #{i} else smoothGroups[sg][i] = true
      )
      if smoothGroups.count == 2 then
      (
         local allfaces;
         local numCaps = 1;
         local theCap;
         local elem1, elem2, elem3;
         local verts1, verts2, verts3;
         local numverts1, numverts2, numverts3;
         
         meshop.detachFaces o smoothGroups[1]
         meshop.detachFaces o smoothGroups[2]
         allfaces = #{1..o.numfaces}
         elem1 = meshop.getElementsUsingFace o 1  -- this method can be improved by PolyTools algorithm
         allfaces -= elem1
         verts1 = (meshop.getVertsUsingFace o elem1) as array; numverts1 = verts1.count 
         elem2 = meshop.getElementsUsingFace o (allfaces as array)[1]
         allfaces -= elem2
         if allfaces.numberset == 0 then      --   Only one cap
         (
            verts2 = (meshop.getVertsUsingFace o elem2) as array; numverts2 = verts2.count
         )
         else
         (
            verts2 = (meshop.getVertsUsingFace o elem2) as array; numverts2 = verts2.count
            
            elem3 = meshop.getElementsUsingFace o (allfaces as array)[1]
            allfaces -= elem3
            if allfaces.numberset != 0 then      --   More elements than expected
            (
               free o
               return false
            )
            else
            (
               numCaps = 2
            )
            verts3 = (meshop.getVertsUsingFace o elem3) as array; numverts3 = verts3.count
         )
         
         if numCaps == 1 then
         (
            theCap = if numverts1 < numverts2 then verts1 else verts2
         )
         else
         (
            myCase = if numverts1 == numverts2 then 1 else if numverts1 == numverts3 then 2 else if numverts2 == numverts3 then 3 else 4
            case myCase of
            (
               1: theCap = verts1
               2: theCap = verts1
               3: theCap = verts2
               4: (free o; return false)   --   Both caps don't have the same number of verts
            )
         )
            
         -- Let's take 3 points of the selected cap (first, middle and last) to find the circle center and Radius
         numV = theCap.count
         p1 =  getVert o theCap[1]; p2 =  getVert o theCap[numV/2]; p3 =  getVert o theCap[numV]
         center = FindCenterOfCircle p1 p2 p3
         radius = distance center p1
         
         numErrors = 0   -- Allow 1 error when checking against cap center if it exists
         for v = 2 to numV-1 do
         (
            if (dd = abs((distance (getVert o theCap[v]) center) - radius)) > 0.01 do numErrors += 1 
            if numErrors > 1 do (free o; return false)   --   More than one distance error. Not a circle
         )
         
         free o
         return true   
      )
      else
      (
         free o
         return false
      )
   )
   
   -- Test:
   for o in (objects as array) do
   (
      itIs = isCylinder o
      if itIs then
      (
         format "Object % is a Cylinder
" o.name
      )
   )
   


Upps! There may be an error if I take the center point of the cap to build the circle…

Deleted… don’t like the solution to avoid using cap center.

Here’s a new version that avoids the mentioned problem:


   
   fn FindCenterOfCircle p1 p2 p3 =
   (
      local a = p3 - p2
      local b = p1 - p3
      local c = p2 - p1
      local u = (dot a a) * (dot c b)
      local v = (dot b b) * (dot c a)
      local w = (dot c c) * (dot b a)
      (u*p1 + v*p2 + w*p3) / (u + v + w)
   )

   
   fn isCylinder theNode =
   (
      if not isKindOf theNode GeometryClass  or isKindOf theNode targetObject do (return false)
      o = snapshotAsMesh theNode
      meshop.autoSmooth o #{1..o.numfaces} 89.9
      smoothGroups = #()
      for i=1 to o.numfaces do
      (
         sg = getFaceSmoothGroup o i
         if smoothGroups[sg] == undefined then smoothGroups[sg] = #{i} else smoothGroups[sg][i] = true
      )
      if smoothGroups.count == 2 then
      (
         local allfaces;
         local numCaps = 1;
         local theCap;
         local elem1, elem2, elem3;
         local verts1, verts2, verts3;
         local numverts1, numverts2, numverts3;
         
         meshop.detachFaces o smoothGroups[1]
         meshop.detachFaces o smoothGroups[2]
         allfaces = #{1..o.numfaces}
         elem1 = meshop.getElementsUsingFace o 1  -- this method can be improved by PolyTools algorithm
         allfaces -= elem1
         verts1 = (meshop.getVertsUsingFace o elem1); numverts1 = verts1.numberSet 
         elem2 = meshop.getElementsUsingFace o (allfaces as array)[1]
         allfaces -= elem2
         if allfaces.numberset == 0 then      --   Only one cap
         (
            verts2 = (meshop.getVertsUsingFace o elem2); numverts2 = verts2.numberSet
         )
         else
         (
            verts2 = (meshop.getVertsUsingFace o elem2); numverts2 = verts2.numberSet
            
            elem3 = meshop.getElementsUsingFace o (allfaces as array)[1]
            allfaces -= elem3
            if allfaces.numberset != 0 then      --   More elements than expected
            (
               free o
               return false
            )
            else
            (
               numCaps = 2
            )
            verts3 = (meshop.getVertsUsingFace o elem3); numverts3 = verts3.numberSet
         )
         
         if numCaps == 1 then
         (
            theCap = if numverts1 < numverts2 then verts1 else verts2
         )
         else
         (
            myCase = if numverts1 == numverts2 then 1 else if numverts1 == numverts3 then 2 else if numverts2 == numverts3 then 3 else 4
            case myCase of
            (
               1: theCap = verts1
               2: theCap = verts1
               3: theCap = verts2
               4: (free o; return false)   --   Both caps don't have the same number of verts
            )
         )
         
         allOpenEdges = meshop.getOpenEdges o
         openVerts = meshop.getVertsUsingEdge o allOpenEdges
         theCap *= openVerts
         theCap = theCap as array
         -- Let's take 3 points of the selected cap (first, middle and last) to find the circle center and Radius
         numV = theCap.count
         p1 =  getVert o theCap[1]; p2 =  getVert o theCap[numV/2]; p3 =  getVert o theCap[numV]
         center = FindCenterOfCircle p1 p2 p3
         radius = distance center p1         
         
         for v = 2 to numV-1 do
         (
            if (dd = abs((distance (getVert o theCap[v]) center) - radius)) > 0.01 do (free o; return false)
         )
         
         free o
         return true   
      )
      else
      (
         free o
         return false
      )
   )
   
   -- Test:
   for o in (objects as array) do
   (
      itIs = isCylinder o
      if itIs then
      (
         format "Object % is a Cylinder
" o.name
      )
   )
   


NOTE: I think all our solutions don’t work with a truncated cone.

As every algorithm based on threshold values, it is very susceptible to them, but this seems to trap small changes in geometry.


(
    
    fn GetVerticesCenter mesh verts =
    (
        cx = cy = cz = 0.0d0
        
        for j in verts do
        (
            vert = getvert mesh j
            cx += vert.x
            cy += vert.y
            cz += vert.z
        )
        
        cx /= verts.numberset
        cy /= verts.numberset
        cz /= verts.numberset
        
        return [cx, cy, cz]
    )
    
    fn IsObjectCylinder node threshold1:0.0001 threshold2:0.05 =
    (
        tmesh = snapshotasmesh node
        faces = tmesh.faces as bitarray
        
        faces1 = #{}
        faces2 = #{}
        
        meshop.autosmooth tmesh faces 75.0

        maxSG = 0
        
        for j in faces do
        (
            sg = getfacesmoothgroup tmesh j
            faces1[j] = sg == 1
            maxSG = amax maxSG sg
        )
        
        if maxSG != 2 do return false
        if (faces-faces1).numberset == 0 do return false
        
        faces2 = faces - faces1
        
        for j in faces1 do
        (
            n1 = getfacenormal tmesh j
            for k in faces2 do
            (
                n2 = getfacenormal tmesh k
                if abs (dot n1 n2) > threshold1 do return false
            )
        )
        
        meshop.detachfaces tmesh faces1
        
        edges = (meshop.getedgesusingface tmesh faces1) * meshop.getopenedges tmesh
        verts = meshop.getvertsusingedge tmesh edges
        
        center = GetVerticesCenter tmesh verts
        
        minDist =  1e9
        maxDist = -1e9
        
        for j in verts do
        (
            d = distance center (getvert tmesh j)
            minDist = amin minDist d
            maxDist = amax maxDist d
        )
        
        return (maxDist-minDist < threshold2)
    )
    
    
    try destroydialog ::RO_TEST catch()
    rollout RO_TEST "" width:192 height:120
    (
        button btn1 "Objects 1" pos:[ 8, 8] width:88 height:32
        button btn2 "Objects 2" pos:[96, 8] width:88 height:32
        button btn3 "Objects 3" pos:[ 8,40] width:88 height:32
        button btn4 "Objects 4" pos:[96,40] width:88 height:32
        button btn5 "Run Test"  pos:[ 8,80] width:176 height:32
        
        local nodes
        
        fn BuildObjects =
        (
            gc()
            delete objects
            
            nodes = for j = 0 to 7 collect
            (
                cylinder heightsegs:10 capsegs:1 sides:(j+3) height:100 radius:10 pos:[j*30,0,0] wirecolor:white
            )
            
            append nodes (cylinder heightsegs:10 capsegs:1 sides:64 height:100 radius:10 pos:[240,0,0] wirecolor:white)
            append nodes (cone heightsegs:10 capsegs:1 sides:10 height:100 radius1:10 radius2:9.98 pos:[270,0,0] wirecolor:white)
            
            converttomesh nodes
        )
        
        fn AddModifiers =
        (
            addmodifier nodes[3] (bend bendangle:0.2)
            addmodifier nodes[4] (noisemodifier strength:[0.05,0.05,0.05])
            addmodifier nodes[5] (push push_value:0.01)
            addmodifier nodes[6] (skew amount:1.0)
            addmodifier nodes[7] (taper amount:0.01)
            addmodifier nodes[8] (twist angle:0.5)
            addmodifier nodes[9] (melt melt_amount:0.5)
        )
        
        fn DeleteCap =
        (
            for j in nodes do deletevert j (j.numverts)
        )
        
        on btn1 pressed do
        (
            BuildObjects()
        )
        
        on btn2 pressed do
        (
            BuildObjects()
            AddModifiers()
        )
        
        on btn3 pressed do
        (
            BuildObjects()
            DeleteCap()
        )
        
        on btn4 pressed do
        (
            BuildObjects()
            DeleteCap()
            AddModifiers()
        )
        
        on btn5 pressed do
        (
            for j in objects do
            (
                isCylinder = IsObjectCylinder j threshold1:0.0001 threshold2:0.05
                if isCylinder then j.wirecolor = green else j.wirecolor = red
                
                format "% -> %
" j.name isCylinder
            )
        )
        
    )
    createdialog RO_TEST
    
)

I’ve waked up this morning with exact and simple formula of the cylinder from the topology view point – all edges are either parallel or perpendicular to the central axis, and all vertices of perpendicular faces are at the same distance from the axis.
another rule:
in every face one edge is parallel and one is perpendicular, or all edges are perpendicular to the central axis.

Here’s a new version based in PolyTools’s one, that avoids checking perpendicularity for all cap faces against all other faces by:

  • checking coplanarity of each cap
  • checking parallelism of both caps (if there are two caps)
    If true, we can check all the other faces just against one cap face, with a really big performance improvement.
    Moreover, with this method we can:
  • Search for the center just of one of the caps
  • Check the radius shape just for one cap. (we can do that as we already have checked coplanarity, parallelism and perpendicularity).

(
   fn GetVerticesCenter p1 p2 p3 =
   (
      local a = p3 - p2
      local b = p1 - p3
      local c = p2 - p1
      local u = (dot a a) * (dot c b)
      local v = (dot b b) * (dot c a)
      local w = (dot c c) * (dot b a)
      (u*p1 + v*p2 + w*p3) / (u + v + w)
   )
   
   fn areCoplanar mesh verts threshold1:0.0001 =
   (
      numV = verts.numberSet
      if numV < 4 do return true
      vertCoords = for v in verts collect getVert mesh v
      p0 = vertCoords[1]
      dir = cross (vertCoords[2] - p0) (vertCoords[3] - p0)
      for i = 4 to numV do
      (
         dotprod = abs(dot dir (vertCoords[i] - p0))
         if dotprod > threshold1 do (return false)      --   format "Exit Not Coplanar: %
" dotprod; 
      )
      true
   )
   
    fn IsObjectCylinder node threshold1:0.0001 threshold2:0.01 =
    (
      local faces, faces1, faces2, faces3
      
        tmesh = snapshotasmesh node
        faces = tmesh.faces as bitarray
      -- Create smooth groups for faces
        meshop.autosmooth tmesh faces 75.0
      smoothGroups = #()
        for j in faces do
        (
            sg = getfacesmoothgroup tmesh j
         if smoothGroups[sg] == undefined then smoothGroups[sg] = #{j} else smoothGroups[sg][j] = true
      )
      if smoothGroups.count != 2 do (free tmesh; return false)      --   print "Exit not 2 smooth groups"; 
   
      
      -- Separate cap faces from others mesh faces
      local numCaps = 1;
      local theCap1, theCap2;
      local elem1, elem2, elem3;
      local verts1, verts2, verts3;
      local numverts1, numverts2, numverts3;
      
      meshop.detachFaces tmesh smoothGroups[1]
   
      elem1 = meshop.getElementsUsingFace tmesh 1  
      faces -= elem1
      verts1 = (meshop.getVertsUsingFace tmesh elem1); numverts1 = verts1.numberSet 
   
      elem2 = meshop.getElementsUsingFace tmesh (faces as array)[1]
      faces -= elem2
      if faces.numberset == 0 then      --   Only one cap
      (
         verts2 = (meshop.getVertsUsingFace tmesh elem2); numverts2 = verts2.numberSet
      )
      else
      (
         verts2 = (meshop.getVertsUsingFace tmesh elem2); numverts2 = verts2.numberSet
         
         elem3 = meshop.getElementsUsingFace tmesh (faces as array)[1]
         faces -= elem3
         if faces.numberset != 0 then      
         (
            free tmesh
            return false                  --   print "Exit More elements than expected"; 
         )
         else
         (
            numCaps = 2
         )
         verts3 = (meshop.getVertsUsingFace tmesh elem3); numverts3 = verts3.numberSet
      )
   
      
      -- store cap faces and cap verts in working variables
      if numCaps == 1 then
      (
         theCap1 = if numverts1 < numverts2 then verts1 else verts2
         if numverts1 < numverts2 then
         (
            theCap1 = verts1; faces1 = elem1; faces3 = elem2
         )
         else
         (
            theCap1 = verts2; faces1 = elem2; faces3 = elem1
         )
      )
      else
      (
         myCase = if numverts1 == numverts2 then 1 else if numverts1 == numverts3 then 2 else if numverts2 == numverts3 then 3 else 4
         case myCase of
         (
            1: (theCap1 = verts1; faces1 = elem1; theCap2 = verts2; faces2 = elem2; faces3 = elem3)
            2: (theCap1 = verts1; faces1 = elem1; theCap2 = verts3; faces2 = elem3; faces3 = elem2)
            3: (theCap1 = verts2; faces1 = elem2; theCap2 = verts3; faces2 = elem3; faces3 = elem1)
            4: (free tmesh; return false)         --   print "Exit Both caps don't have the same number of verts";
         )
      )
      
      -- Check for caps coplanarity and parallelism
      local first = (faces1 as array)[1]
      local second   
      local third = (faces3 as array)[1]
      
      coplanar = areCoplanar tmesh theCap1
      if not coplanar do (free tmesh; return false)      --   print "Exit Cap1 not coplanar"; 
      if faces2 != undefined do
      (
         coplanar = areCoplanar tmesh theCap2
         if not coplanar do (free tmesh; return false)   --   print "Exit Cap2 not coplanar"; 
         
         second = (faces2 as array)[1]
         parallel = length (cross (getfacenormal tmesh first) (getfacenormal tmesh second))
         if parallel > threshold1 do (free tmesh; return false)      --   print "Exit Caps are not parallel"; 
      )
      
      -- Check for perpendicularity between caps and other faces
      n1 = getfacenormal tmesh first 
      for k in faces3 do
      (
         n2 = getfacenormal tmesh k
         if abs (dot n1 n2) > threshold1 do (free tmesh; return false)      --   print "Exit other faces are not perpendicular"; 
      )
      
      -- Check for circular shape of caps
        edges = (meshop.getedgesusingface tmesh faces1) * meshop.getopenedges tmesh
        verts = meshop.getvertsusingedge tmesh edges
      vertsArray = verts as array
      numV = vertsArray.count
      p1 =  getVert tmesh vertsArray[1]; p2 =  getVert tmesh vertsArray[numV/3+1]; p3 =  getVert tmesh vertsArray[numV*2/3+1]
      center = GetVerticesCenter p1 p2 p3
        minDist =  1e9
        maxDist = -1e9
        for j in verts do
        (
            d = distance center (getvert tmesh j)
            minDist = amin minDist d
            maxDist = amax maxDist d
        )
      free tmesh; 
      --format "minDist= %; maxDist= %
" minDist maxDist
        return (maxDist-minDist < threshold2)
    )
   
   
    try destroydialog ::RO_TEST catch()
    rollout RO_TEST "" width:192 height:120
    (
        button btn1 "Objects 1" pos:[ 8, 8] width:88 height:32
        button btn2 "Objects 2" pos:[96, 8] width:88 height:32
        button btn3 "Objects 3" pos:[ 8,40] width:88 height:32
        button btn4 "Objects 4" pos:[96,40] width:88 height:32
        button btn5 "Run Test"  pos:[ 8,80] width:176 height:32
        local nodes
        fn BuildObjects =
        (
            gc()
            delete objects
            nodes = for j = 0 to 7 collect
            (
                cylinder heightsegs:10 capsegs:1 sides:(j+3) height:100 radius:10 pos:[j*30,0,0] wirecolor:white
            )
            append nodes (cylinder heightsegs:10 capsegs:1 sides:64 height:100 radius:10 pos:[240,0,0] wirecolor:white)
            append nodes (cone heightsegs:10 capsegs:1 sides:10 height:100 radius1:10 radius2:9.98 pos:[270,0,0] wirecolor:white)
            converttomesh nodes
        )
        fn AddModifiers =
        (
            addmodifier nodes[3] (bend bendangle:0.2)
            addmodifier nodes[4] (noisemodifier strength:[0.05,0.05,0.05])
            addmodifier nodes[5] (push push_value:0.01)
            addmodifier nodes[6] (skew amount:1.0)
            addmodifier nodes[7] (taper amount:0.01)
            addmodifier nodes[8] (twist angle:0.5)
            addmodifier nodes[9] (melt melt_amount:0.5)
        )
        fn DeleteCap =
        (
            for j in nodes do deletevert j (j.numverts)
        )
        on btn1 pressed do
        (
            BuildObjects()
        )
        on btn2 pressed do
        (
            BuildObjects()
            AddModifiers()
        )
        on btn3 pressed do
        (
            BuildObjects()
            DeleteCap()
        )
        on btn4 pressed do
        (
            BuildObjects()
            DeleteCap()
            AddModifiers()
        )
        on btn5 pressed do
        (
         t1=timestamp()
            for j in objects do
            (
                isCylinder = IsObjectCylinder j threshold1:0.0001 threshold2:0.01
                if isCylinder then j.wirecolor = green else j.wirecolor = red
                format "% -> %
" j.name isCylinder
            )
         t2=timestamp()
         format "Time: %ms
" (t2-t1)
        )
    )
    createdialog RO_TEST
)
   


Credits:

  • @miauu for the interesting post. What will he do with that!!??
  • @Serejah for his brilliant idea of using smooth groups (how did you came up to this??!!)
  • @PolyTools3D, an algorithms crack
  • @Swordslayer for his circumcenter function
  • @DenisT and @vusta for their contributions
  • @aaandres, just a copy-paster
Page 2 / 4