Search Unity

How do I get the cross sectional surface area of my geometry mesh in meters?

Discussion in 'General Graphics' started by AlanMattano, Mar 16, 2017.

  1. AlanMattano

    AlanMattano

    Joined:
    Aug 22, 2013
    Posts:
    1,501
    Last edited: Mar 21, 2017
  2. AlanMattano

    AlanMattano

    Joined:
    Aug 22, 2013
    Posts:
    1,501
    Probably a difficult question for an expert.
     
  3. MatthewFK

    MatthewFK

    Joined:
    Jul 5, 2021
    Posts:
    1
    I am using this code for an airplane. It assumes that your models are convex. If you want concave models to work, you will have to break up your model (and modify my algorithm to stop including all of the child objects when finding vertices).

    The logic goes like this:

    Get an array of Vector3s containing all of the vertices (that are different). So that we do not have to find this every time, I store it in a dictionary.

    It finds two perpendicular vectors (perp1 and perp2) that are perpendicular to the normal vector. These will be used as basis vectors.

    It squishes all of the vertices on the plane perpendicular to normal using the Dot method to find these points in terms of the basis vectors in the plane.

    It uses a hull mapping algorithm from https://gist.github.com/dLopreiato/7fd142d0b9728518552188794b8a750c to find the hull of the shape.

    It finds the area of the hull and returns it.

    Code (CSharp):
    1. public static Dictionary<GameObject, Vector3[]> verticesDict = new Dictionary<GameObject, Vector3[]>();
    2.  
    3.     public static float CrossSectionArea(GameObject g, Vector3 normal)
    4.     {
    5.         if(normal.x == 0f && normal.y == 0f && normal.z == 0f)
    6.         {
    7.             return 0f;
    8.         }
    9.  
    10.         Vector3[] vertices;
    11.        
    12.         PutVerticesInDict(g, out vertices);
    13.        
    14.  
    15.         normal = normal.normalized;
    16.  
    17.         Vector3 perp1 = Perp(normal).normalized;
    18.         Vector3 perp2 = Vector3.Cross(normal, perp1);
    19.  
    20.  
    21.         List<Vector2> inThePlane = new List<Vector2>(vertices.Length);
    22.        
    23.  
    24.         for (int i = 0; i < vertices.Length; i++)
    25.         {
    26.            
    27.             inThePlane.Add(new Vector2(Vector3.Dot(perp1, vertices[i]), Vector3.Dot(perp2, vertices[i])));
    28.         }
    29.         return HullArea(ComputeConvexHull(inThePlane));
    30.     }
    31.     public static void PutVerticesInDict(GameObject g, out Vector3[] vertices)
    32.     {
    33.         if (!verticesDict.TryGetValue(g, out vertices))
    34.         {
    35.             MeshFilter[] meshes = g.GetComponentsInChildren<MeshFilter>();
    36.  
    37.  
    38.             HashSet<Vector3> svertices = new HashSet<Vector3>();
    39.             foreach (MeshFilter m in meshes)
    40.             {
    41.  
    42.                 foreach (Vector3 v3 in m.mesh.vertices)
    43.                 {
    44.  
    45.                     svertices.Add(m.transform.TransformPoint(v3));
    46.                 }
    47.             }
    48.  
    49.             vertices = new Vector3[svertices.Count];
    50.             svertices.CopyTo(vertices);
    51.             verticesDict.Add(g, vertices);
    52.         }
    53.     }
    54.     public static Vector3 Perp(Vector3 v3)
    55.     {
    56.         return v3.z < v3.x ? new Vector3(v3.y, -v3.x, 0) : new Vector3(0, -v3.z, v3.y);
    57.     }
    58.  
    59.  
    60.     public static float HullArea(IList<Vector2> hull)
    61.     {
    62.        
    63.         float sum = 0f;
    64.         for (int i = 1; i < hull.Count - 1; i++)
    65.         {
    66.             sum += TriangleArea(hull[0], hull[i], hull[i + 1]);
    67.         }
    68.         return sum;
    69.     }
    70.     public static float TriangleArea(Vector2 v0, Vector2 v1, Vector2 v2)
    71.     {
    72.  
    73.         return Mathf.Abs((v1.x - v0.x) * (v2.y - v0.y) - (v2.x - v0.x) * (v1.y - v0.y)) * .5f;
    74.     }
    75.     //From https://gist.github.com/dLopreiato/7fd142d0b9728518552188794b8a750c
    76.     public static IList<Vector2> ComputeConvexHull(List<Vector2> points, bool sortInPlace = false)
    77.     {
    78.         if (!sortInPlace)
    79.             points = new List<Vector2>(points);
    80.         points.Sort((a, b) =>
    81.             a.x == b.x ? a.y.CompareTo(b.y) : (a.x > b.x ? 1 : -1));
    82.  
    83.         // Importantly, DList provides O(1) insertion at beginning and end
    84.         CircularList<Vector2> hull = new CircularList<Vector2>();
    85.         int L = 0, U = 0; // size of lower and upper hulls
    86.  
    87.         // Builds a hull such that the output polygon starts at the leftmost Vector2.
    88.         for (int i = points.Count - 1; i >= 0; i--)
    89.         {
    90.             Vector2 p = points[i], p1;
    91.  
    92.             // build lower hull (at end of output list)
    93.             while (L >= 2 && ((p1 = hull.Last) - (hull[hull.Count - 2])).Cross(p - p1) >= 0)
    94.             {
    95.                 hull.PopLast();
    96.                 L--;
    97.             }
    98.             hull.PushLast(p);
    99.             L++;
    100.  
    101.             // build upper hull (at beginning of output list)
    102.             while (U >= 2 && ((p1 = hull.First) - (hull[1])).Cross(p - p1) <= 0)
    103.             {
    104.                 hull.PopFirst();
    105.                 U--;
    106.             }
    107.             if (U != 0) // when U=0, share the Vector2 added above
    108.                 hull.PushFirst(p);
    109.             U++;
    110.             Debug.Assert(U + L == hull.Count + 1);
    111.         }
    112.         hull.PopLast();
    113.         return hull;
    114.     }
    115.  
    116.     private static Vector2 Sub(this Vector2 a, Vector2 b)
    117.     {
    118.         return a - b;
    119.     }
    120.  
    121.     private static float Cross(this Vector2 a, Vector2 b)
    122.     {
    123.         return a.x * b.y - a.y * b.x;
    124.     }
    125.  
    126.     private class CircularList<T> : List<T>
    127.     {
    128.         public T Last
    129.         {
    130.             get
    131.             {
    132.                 return this[this.Count - 1];
    133.             }
    134.             set
    135.             {
    136.                 this[this.Count - 1] = value;
    137.             }
    138.         }
    139.  
    140.         public T First
    141.         {
    142.             get
    143.             {
    144.                 return this[0];
    145.             }
    146.             set
    147.             {
    148.                 this[0] = value;
    149.             }
    150.         }
    151.  
    152.         public void PushLast(T obj)
    153.         {
    154.             this.Add(obj);
    155.         }
    156.  
    157.         public T PopLast()
    158.         {
    159.             T retVal = this[this.Count - 1];
    160.             this.RemoveAt(this.Count - 1);
    161.             return retVal;
    162.         }
    163.  
    164.         public void PushFirst(T obj)
    165.         {
    166.             this.Insert(0, obj);
    167.         }
    168.  
    169.         public T PopFirst()
    170.         {
    171.             T retVal = this[0];
    172.             this.RemoveAt(0);
    173.             return retVal;
    174.         }
    175.     }
     
    Last edited: Oct 13, 2021
  4. AlanMattano

    AlanMattano

    Joined:
    Aug 22, 2013
    Posts:
    1,501
    Hello Matthew,
    Your first answer is amazing. And by the way, very well come to this forum!

    I read your code and I did not implement it yet.
    But since I'm relatively new to C#, I did not understand
    Code (CSharp):
    1. Cross(p - p1)
    at lines 93 and 102 because Cross needs two inputs.

    I tried rewriting but I came across this line
    Code (CSharp):
    1. ((p1 = hull.Last) - (hull[hull.Count - 2])).Cross(p - p1)
    Can it be then instead of a - is a ","? for exmple
    Code (CSharp):
    1. // build lower hull (at end of output list)
    2.             while ( L >= 2 )
    3.             {
    4.                 p1 = hull.Last - (hull[hull.Count - 2]);
    5.  
    6.                 if ( Cross(p , p1) >= 0)
    7.                 {
    8.                     hull.PopLast();
    9.                     L--;
    10.                 }
    11.             }
     
    Last edited: Oct 25, 2021
  5. colinalaws

    colinalaws

    Joined:
    Sep 7, 2019
    Posts:
    4
    I have an idea that would work pretty decently as an approximation possibly. Maybe not something for seriously realistic game. If the cross sectional area is really just a 2D area looking from the front of the object, then why not just place some transforms into the world around each component or mesh, and then at runtime compute the area of the shape you generated from 2D vectors.

    This of course relies on (my limited) understanding of this formula. If in fact the surface area of the front of the mesh is needed, then certainly you would need to iterate over the area and compute this. I would recommend performing this once for your parts, store it as a const, and use it in an AirResistanceModifier (or something named like that) MonoBehaviour that updates the drag property of the rigidbody component.

    I'm needing to implement this right now, and I figured, hey why not just place down some vectors to define the general shape of the object. If you're looking for simulation, this is not what you want of course. Anyways, thought I'd mention my idea, since I'm working on the same thing.
     
  6. koirat

    koirat

    Joined:
    Jul 7, 2012
    Posts:
    2,073
    Since we are talking about drag we do not want a cross section, what we need is a projection of a mesh in direction of fluid movement.

    An approximate solution could be obtained by rendering mesh on a texture and than counting number of pixels. Can be calculated by montecarlo method for faster computations.
     
    colinalaws likes this.
  7. colinalaws

    colinalaws

    Joined:
    Sep 7, 2019
    Posts:
    4
    Gotcha! Thank you so much for clarifying this confusing bit for me. Examples of cross sectional area always resulted in images of 2D planes intersecting 3D objects. Your solution makes perfect sense to me, and I can definitely implement that. Thank you very much for helping me find a simple way to easily compute this.
     
  8. colinalaws

    colinalaws

    Joined:
    Sep 7, 2019
    Posts:
    4
    I have two ideas here after studying the Monte Carlo method.

    1) Place a camera in the scene in orthographic mode which only renders the flying craft. Then, the camera would be programmed to move in the direction that the fluid is moving, and use the Monte Carlo method to compute the cross sectional area for a given frame that was captured. I know I'd need some way to scale the measurement coming from the camera, essentially I'd need to know the number of units my viewport is measuring.

    2) This is just me taking a stab in the dark, as I am not well versed in this subject. But I wonder if I'd be possible to create a heuristic model to determine drag given the direction the fluid would be moving over the craft. I was watching a video about Monte Carlo tree search, and thought maybe it could quickly calculate the drag. The applications where for long task decision making, but I feel like the application could be much broader.
     
  9. SebLazyWizard

    SebLazyWizard

    Joined:
    Jun 15, 2018
    Posts:
    234
    That's an interesting topic and the idea to count pixels of a mesh rendered to a texture is a great method to get the frontal area.
    But to properly simulate drag you also have to include the drag coefficient of the mesh;
    https://en.wikipedia.org/wiki/Drag_coefficient

    Maybe you can render the camera's depth to a texture and analyse the meshes shape to somehow approximate the drag coefficient.
    Like building up a curvature map from it to get convex/concave info, or even create a flowmap based on that.
     
    Last edited: Jan 29, 2022