Attribute VB_Name = "Geometry_Cell"

Option Explicit

'Tolerances
Public Const Geometry_WeldTolerance As Double = 0.00001 'How close to a position is the same position.
Public Const Geometry_OnPlaneTolerance As Double = Geometry_WeldTolerance * 10 'How close to a surface plane is on the plane.

'Limits
Public Const Geometry_MaxPoints As Long = 64
Public Const Geometry_MaxSufraces As Long = 32
Public Const Geometry_MaxIndices As Long = 30

'Surface types
Public Const Geometry_SurfaceType_Default As Byte = 0
Public Const Geometry_SurfaceType_Infinity As Byte = 1

'Surfaces
Public Type Surface
    'いい Place shallow material index here (Actually displayed on the surface)
    'いい Set visibility here to create windows and save area on lightmaps
    'いい Store a direction for texture coordinate generation so that wood can look correct no matter what direction
    'Primary
    SurfacePlane As FractionPlane
    'Secondary
    SurfaceType As Byte 'Geometry_SurfaceType_Default or Geometry_SurfaceType_Infinity
    NumberOfIndices As Byte
    Indices(0 To Geometry_MaxIndices - 1) As Byte
End Type

'Cells
Public Type CellSettings
    MaterialIndex As Long 'いい Use to select a material
    Selected As Boolean
    IsBroken As Boolean
End Type
'Invariant:
'   Each cell must have at least one surface
'   Each cell has a volume greater than zero
'   Each cell is convex
Public Type Cell
    'Stored
    Surfaces(0 To Geometry_MaxSufraces - 1) As Surface
    NumberOfSurfaces As Byte
    'Temporary
    Positions(0 To Geometry_MaxPoints - 1) As Vector3
    NumberOfPositions As Byte
    Bound As Box 'いい Use for broad phases
    'Inherited data
    Settings As CellSettings
End Type
Public DefaultCell As Cell
Public Type CellIntersectionResult
    PlaneIntersection As PlaneIntersectionResult
    CellIndex As Long
    SurfaceIndex As Long
End Type
Public Type CellToPlaneResult
    Fronts As Integer
    Centers As Integer
    Backs As Integer
End Type

Public Function MakeCellSettings(MaterialIndex As Long, Optional Selected As Boolean = False, Optional IsBroken As Boolean = False) As CellSettings
    MakeCellSettings.MaterialIndex = MaterialIndex
    MakeCellSettings.Selected = Selected
    MakeCellSettings.IsBroken = IsBroken
End Function

Public Function ValidateCell(C As Cell, ByRef Errors As String) As Boolean
    Dim SurfaceIndex As Long
    Dim PointCount As Long
    Dim PointIndex As Long
    ValidateCell = True
    For SurfaceIndex = 0 To C.NumberOfSurfaces - 1
        PointCount = C.Surfaces(SurfaceIndex).NumberOfIndices
        If PointCount < 3 Then
            Errors = Errors & "Not enough points to define a normal in surface #" & SurfaceIndex & vbNewLine
            ValidateCell = False
        End If
        For PointIndex = 0 To PointCount - 1
            Dim VertexIndex As Long
            Dim NextVertexIndex As Long
            Dim NextPointIndex As Long
            Dim Point As Vector3
            Dim SurfacePlane As Plane
            Dim Dist As Double
            VertexIndex = C.Surfaces(SurfaceIndex).Indices(PointIndex)
            NextPointIndex = (PointIndex + 1) Mod PointCount
            NextVertexIndex = C.Surfaces(SurfaceIndex).Indices(NextPointIndex)
            If VertexIndex = NextVertexIndex Then
                Errors = Errors & "Multiple points " & PointIndex & " and " & NextPointIndex & " in surface #" & SurfaceIndex & " have the same index #" & VertexIndex & vbNewLine
                ValidateCell = False
            End If
            Point = C.Positions(VertexIndex)
            SurfacePlane = ApproximateFractionPlane(C.Surfaces(SurfaceIndex).SurfacePlane)
            Dist = PointToSurfaceDistance(Point, SurfacePlane)
            If Abs(Dist) > Geometry_OnPlaneTolerance Then
                Errors = Errors & "Point not at plane (" & Dist & ") in surface #" & SurfaceIndex & ", vertex #" & PointIndex & vbNewLine
                ValidateCell = False
            End If
        Next PointIndex
    Next SurfaceIndex
End Function

Public Sub UpdateCellBound(C As Cell)
    If C.NumberOfPositions > 0 Then
        Dim Bound As Box
        Dim P As Long
        Bound = MakeBoxFromTwoPoints(C.Positions(0), C.Positions(0))
        For P = 1 To C.NumberOfPositions - 1
            Bound = ExtendBoxUsingPoint(Bound, C.Positions(P))
        Next P
        C.Bound = Bound
    End If
End Sub

Public Function GenerateCellFromPlanes(C As Cell) As Cell
    Dim CurrentCell As Cell
    CurrentCell = MakeCell_Box(MakeFractionBox_Big, C.Settings, True)
    Dim SurfaceIndex As Long
    For SurfaceIndex = 0 To C.NumberOfSurfaces - 1
        If C.Surfaces(SurfaceIndex).SurfaceType <> Geometry_SurfaceType_Infinity And CellIntersectsPlane(CurrentCell, C.Surfaces(SurfaceIndex).SurfacePlane) Then
            CurrentCell = GetFragmentOfCell(CurrentCell, FlipFractionPlane(C.Surfaces(SurfaceIndex).SurfacePlane))
        End If
    Next SurfaceIndex
    GenerateCellFromPlanes = CurrentCell
End Function

Public Function GetCornerDirection(C As Cell, PointIndex As Long) As Vector3
    Dim Planes(0 To Geometry_MaxSufraces - 1) As Plane
    Dim NumberOfPlanes As Long
    Dim SurfaceIndex As Long
    Dim P As Long
    Dim I As Long
    For SurfaceIndex = 0 To C.NumberOfSurfaces - 1
        For P = 0 To C.Surfaces(SurfaceIndex).NumberOfIndices - 1
            If C.Surfaces(SurfaceIndex).Indices(P) = PointIndex Then
                Planes(NumberOfPlanes) = ApproximateFractionPlane(C.Surfaces(SurfaceIndex).SurfacePlane)
                NumberOfPlanes = NumberOfPlanes + 1
                Exit For
            End If
        Next P
    Next SurfaceIndex
    If NumberOfPlanes < 3 Then
        GetCornerDirection = MakeVector3(0, 0, 0)
        Exit Function
    End If
    Dim StartLocation As Vector3
    Dim CurrentLocation As Vector3
    StartLocation = C.Positions(PointIndex)
    CurrentLocation = StartLocation
    For I = 1 To 32
        Dim Error As Double
        Dim Dist As Double
        Dim Correction As Double
        Error = 0
        For P = 0 To NumberOfPlanes - 1
            Dist = 1 - PointToSurfaceDistance(CurrentLocation, Planes(P))
            Correction = Dist * 0.9
            CurrentLocation = AddVector3(CurrentLocation, MulVector3(Planes(P).Normal, Correction))
            Error = Max_Double(Error, Abs(Dist))
        Next P
        If Error < 0.00000001 Then Exit For
    Next I
    GetCornerDirection = ClampSphere3(SubVector3(CurrentLocation, StartLocation), 2.5)
End Function

Public Function MakeCell_Infinite(Settings As CellSettings) As Cell
    MakeCell_Infinite.Settings = Settings
End Function

Public Function MakeCell_Box(B As FractionBox, Settings As CellSettings, IsEnd As Boolean) As Cell
    Dim Indices(0 To 7) As Long
    Indices(0) = AddPosition(MakeCell_Box, ApproximateFractionPoint(MakeFractionPoint_Args(B.Max.X, B.Max.Y, B.Max.Z))) '0 +++
    Indices(1) = AddPosition(MakeCell_Box, ApproximateFractionPoint(MakeFractionPoint_Args(B.Min.X, B.Max.Y, B.Max.Z))) '1 -++
    Indices(2) = AddPosition(MakeCell_Box, ApproximateFractionPoint(MakeFractionPoint_Args(B.Min.X, B.Min.Y, B.Max.Z))) '2 --+
    Indices(3) = AddPosition(MakeCell_Box, ApproximateFractionPoint(MakeFractionPoint_Args(B.Max.X, B.Min.Y, B.Max.Z))) '3 +-+
    Indices(4) = AddPosition(MakeCell_Box, ApproximateFractionPoint(MakeFractionPoint_Args(B.Max.X, B.Min.Y, B.Min.Z))) '4 +--
    Indices(5) = AddPosition(MakeCell_Box, ApproximateFractionPoint(MakeFractionPoint_Args(B.Min.X, B.Min.Y, B.Min.Z))) '5 ---
    Indices(6) = AddPosition(MakeCell_Box, ApproximateFractionPoint(MakeFractionPoint_Args(B.Min.X, B.Max.Y, B.Min.Z))) '6 -+-
    Indices(7) = AddPosition(MakeCell_Box, ApproximateFractionPoint(MakeFractionPoint_Args(B.Max.X, B.Max.Y, B.Min.Z))) '7 ++-
    MakeCell_Box.NumberOfSurfaces = 6
    Dim SurfaceIndex As Long
    SurfaceIndex = 0
    MakeCell_Box.Surfaces(SurfaceIndex).NumberOfIndices = 4 '+Y
    MakeCell_Box.Surfaces(SurfaceIndex).SurfacePlane = MakeFractionPlane(MakeFractionPoint_Args(B.Min.X, B.Max.Y, B.Min.Z), MakeFractionNormal_Args(0, 1, 0))
    MakeCell_Box.Surfaces(SurfaceIndex).Indices(0) = Indices(7)
    MakeCell_Box.Surfaces(SurfaceIndex).Indices(1) = Indices(6)
    MakeCell_Box.Surfaces(SurfaceIndex).Indices(2) = Indices(1)
    MakeCell_Box.Surfaces(SurfaceIndex).Indices(3) = Indices(0)
    SurfaceIndex = 1
    MakeCell_Box.Surfaces(SurfaceIndex).NumberOfIndices = 4 '-Y
    MakeCell_Box.Surfaces(SurfaceIndex).SurfacePlane = MakeFractionPlane(MakeFractionPoint_Args(B.Min.X, B.Min.Y, B.Min.Z), MakeFractionNormal_Args(0, -1, 0))
    MakeCell_Box.Surfaces(SurfaceIndex).Indices(0) = Indices(5)
    MakeCell_Box.Surfaces(SurfaceIndex).Indices(1) = Indices(4)
    MakeCell_Box.Surfaces(SurfaceIndex).Indices(2) = Indices(3)
    MakeCell_Box.Surfaces(SurfaceIndex).Indices(3) = Indices(2)
    SurfaceIndex = 2
    MakeCell_Box.Surfaces(SurfaceIndex).NumberOfIndices = 4 '+X
    MakeCell_Box.Surfaces(SurfaceIndex).SurfacePlane = MakeFractionPlane(MakeFractionPoint_Args(B.Max.X, B.Min.Y, B.Min.Z), MakeFractionNormal_Args(1, 0, 0))
    MakeCell_Box.Surfaces(SurfaceIndex).Indices(0) = Indices(0)
    MakeCell_Box.Surfaces(SurfaceIndex).Indices(1) = Indices(3)
    MakeCell_Box.Surfaces(SurfaceIndex).Indices(2) = Indices(4)
    MakeCell_Box.Surfaces(SurfaceIndex).Indices(3) = Indices(7)
    SurfaceIndex = 3
    MakeCell_Box.Surfaces(SurfaceIndex).NumberOfIndices = 4 '-X
    MakeCell_Box.Surfaces(SurfaceIndex).SurfacePlane = MakeFractionPlane(MakeFractionPoint_Args(B.Min.X, B.Min.Y, B.Min.Z), MakeFractionNormal_Args(-1, 0, 0))
    MakeCell_Box.Surfaces(SurfaceIndex).Indices(0) = Indices(6)
    MakeCell_Box.Surfaces(SurfaceIndex).Indices(1) = Indices(5)
    MakeCell_Box.Surfaces(SurfaceIndex).Indices(2) = Indices(2)
    MakeCell_Box.Surfaces(SurfaceIndex).Indices(3) = Indices(1)
    SurfaceIndex = 4
    MakeCell_Box.Surfaces(SurfaceIndex).NumberOfIndices = 4 '+Z
    MakeCell_Box.Surfaces(SurfaceIndex).SurfacePlane = MakeFractionPlane(MakeFractionPoint_Args(B.Min.X, B.Min.Y, B.Max.Z), MakeFractionNormal_Args(0, 0, 1))
    MakeCell_Box.Surfaces(SurfaceIndex).Indices(0) = Indices(0)
    MakeCell_Box.Surfaces(SurfaceIndex).Indices(1) = Indices(1)
    MakeCell_Box.Surfaces(SurfaceIndex).Indices(2) = Indices(2)
    MakeCell_Box.Surfaces(SurfaceIndex).Indices(3) = Indices(3)
    SurfaceIndex = 5
    MakeCell_Box.Surfaces(SurfaceIndex).NumberOfIndices = 4 '-Z
    MakeCell_Box.Surfaces(SurfaceIndex).SurfacePlane = MakeFractionPlane(MakeFractionPoint_Args(B.Min.X, B.Min.Y, B.Min.Z), MakeFractionNormal_Args(0, 0, -1))
    MakeCell_Box.Surfaces(SurfaceIndex).Indices(0) = Indices(4)
    MakeCell_Box.Surfaces(SurfaceIndex).Indices(1) = Indices(5)
    MakeCell_Box.Surfaces(SurfaceIndex).Indices(2) = Indices(6)
    MakeCell_Box.Surfaces(SurfaceIndex).Indices(3) = Indices(7)
    For SurfaceIndex = 0 To 5
        If IsEnd Then
            MakeCell_Box.Surfaces(SurfaceIndex).SurfaceType = Geometry_SurfaceType_Infinity
        Else
            MakeCell_Box.Surfaces(SurfaceIndex).SurfaceType = Geometry_SurfaceType_Default
        End If
    Next SurfaceIndex
    MakeCell_Box.Settings = Settings
    UpdateCellBound MakeCell_Box
End Function

'Postcondition: Returns the index of the existing position at Position or -1 if not found within Geometry_WeldTolerance
Public Function GetPositionIndex(ByRef TargetCell As Cell, Position As Vector3) As Long
    GetPositionIndex = -1
    Dim Closest As Double
    Closest = Geometry_WeldTolerance
    Dim P As Long
    For P = 0 To TargetCell.NumberOfPositions - 1
        Dim Dist As Single
        Dist = AbsVector3(SubVector3(Position, TargetCell.Positions(P)))
        If Dist < Closest Then
            GetPositionIndex = P
            Closest = Dist
        End If
    Next P
End Function

Public Function AddPosition(ByRef TargetCell As Cell, NewPosition As Vector3) As Long
    Dim ExistingIndex As Long
    ExistingIndex = GetPositionIndex(TargetCell, NewPosition)
    If ExistingIndex > -1 Then
        AddPosition = ExistingIndex
    Else
        TargetCell.Positions(TargetCell.NumberOfPositions) = NewPosition
        AddPosition = TargetCell.NumberOfPositions
        TargetCell.NumberOfPositions = TargetCell.NumberOfPositions + 1
    End If
End Function

Public Function GetBoundingBox(C As Cell) As Box
    Dim SurfaceIndex As Long
    Dim P As Long
    Dim Pos As Vector3
    GetBoundingBox.Min = MakeVector3(1E+30, 1E+30, 1E+30)
    GetBoundingBox.Max = MakeVector3(-1E+30, -1E+30, -1E+30)
    For SurfaceIndex = 0 To C.NumberOfSurfaces - 1
        For P = 0 To C.Surfaces(SurfaceIndex).NumberOfIndices - 1
            Pos = C.Positions(C.Surfaces(SurfaceIndex).Indices(P))
            GetBoundingBox.Min = Min_Vector3(GetBoundingBox.Min, Pos)
            GetBoundingBox.Max = Max_Vector3(GetBoundingBox.Max, Pos)
        Next P
    Next SurfaceIndex
End Function

'Postcondition: True iff C has at least one used point in each corner of the axis aligned bounding box
Public Function IsCellBox(C As Cell, Tolerance As Double) As Boolean
    Dim Bound As Box
    Dim SurfaceIndex As Long
    Dim P As Long
    Dim Pos As Vector3
    Dim Corner(0 To 7) As Boolean
    Bound = GetBoundingBox(C)
    For SurfaceIndex = 0 To C.NumberOfSurfaces - 1
        For P = 0 To C.Surfaces(SurfaceIndex).NumberOfIndices - 1
            Pos = C.Positions(C.Surfaces(SurfaceIndex).Indices(P))
            If AbsVector3(SubVector3(Pos, MakeVector3(Bound.Min.X, Bound.Min.Y, Bound.Min.Z))) < Tolerance Then Corner(0) = True
            If AbsVector3(SubVector3(Pos, MakeVector3(Bound.Max.X, Bound.Min.Y, Bound.Min.Z))) < Tolerance Then Corner(1) = True
            If AbsVector3(SubVector3(Pos, MakeVector3(Bound.Min.X, Bound.Max.Y, Bound.Min.Z))) < Tolerance Then Corner(2) = True
            If AbsVector3(SubVector3(Pos, MakeVector3(Bound.Max.X, Bound.Max.Y, Bound.Min.Z))) < Tolerance Then Corner(3) = True
            If AbsVector3(SubVector3(Pos, MakeVector3(Bound.Min.X, Bound.Min.Y, Bound.Max.Z))) < Tolerance Then Corner(4) = True
            If AbsVector3(SubVector3(Pos, MakeVector3(Bound.Max.X, Bound.Min.Y, Bound.Max.Z))) < Tolerance Then Corner(5) = True
            If AbsVector3(SubVector3(Pos, MakeVector3(Bound.Min.X, Bound.Max.Y, Bound.Max.Z))) < Tolerance Then Corner(6) = True
            If AbsVector3(SubVector3(Pos, MakeVector3(Bound.Max.X, Bound.Max.Y, Bound.Max.Z))) < Tolerance Then Corner(7) = True
        Next P
    Next SurfaceIndex
    IsCellBox = Corner(0) And Corner(1) And Corner(2) And Corner(3) And Corner(4) And Corner(5) And Corner(6) And Corner(7)
End Function

Public Function CellToPlane(C As Cell, CuttingPlane As Plane) As CellToPlaneResult
    Dim P As Long
    Dim Pos As Vector3
    Dim Dist As Double
    Dim HasFront As Boolean
    Dim HasBack As Boolean
    Dim Centers As Integer
    For P = 0 To C.NumberOfPositions - 1
        Pos = C.Positions(P)
        Dist = PointToSurfaceDistance(Pos, CuttingPlane)
        If Dist > Geometry_OnPlaneTolerance Then
            CellToPlane.Fronts = CellToPlane.Fronts + 1
        ElseIf Dist < -Geometry_OnPlaneTolerance Then
            CellToPlane.Backs = CellToPlane.Backs + 1
        Else
            CellToPlane.Centers = CellToPlane.Centers + 1
        End If
    Next P
End Function

Public Function CellIntersectsPlane(C As Cell, CuttingPlane As FractionPlane) As Boolean
    Dim Result As CellToPlaneResult
    Result = CellToPlane(C, ApproximateFractionPlane(CuttingPlane))
    CellIntersectsPlane = (Result.Fronts > 0 And Result.Backs > 0)
End Function

Public Function CellOnPlanePositive(C As Cell, CuttingPlane As FractionPlane) As Boolean
    Dim Result As CellToPlaneResult
    Result = CellToPlane(C, ApproximateFractionPlane(CuttingPlane))
    CellOnPlanePositive = (Result.Fronts > 0 And Result.Backs = 0)
End Function

Public Function CellOnPlaneNegative(C As Cell, CuttingPlane As FractionPlane) As Boolean
    Dim Result As CellToPlaneResult
    Result = CellToPlane(C, ApproximateFractionPlane(CuttingPlane))
    CellOnPlaneNegative = (Result.Fronts = 0 And Result.Backs > 0)
End Function

Public Function CellTouchesPlane(C As Cell, CuttingPlane As FractionPlane) As Boolean
    Dim Result As CellToPlaneResult
    Result = CellToPlane(C, ApproximateFractionPlane(CuttingPlane))
    CellTouchesPlane = (Result.Fronts > 0 And Result.Centers >= 3 And Result.Backs = 0)
End Function

Public Function PointInsideCell(Point As Vector3, C As Cell) As Boolean
    Dim SurfaceIndex As Long
    PointInsideCell = True 'True unless any point is outside
    For SurfaceIndex = 0 To C.NumberOfSurfaces - 1
        If PointToSurfaceDistance(Point, ApproximateFractionPlane(C.Surfaces(SurfaceIndex).SurfacePlane)) > Geometry_OnPlaneTolerance Then
            PointInsideCell = False
            Exit For
        End If
    Next SurfaceIndex
End Function

Public Function LineToCellIntersection(C As Cell, CellIndex As Long, L As Line3) As CellIntersectionResult
    Dim SurfaceIndex As Long
    Dim S2 As Long
    Dim Hit As PlaneIntersectionResult
    LineToCellIntersection.PlaneIntersection.Collided = False
    LineToCellIntersection.CellIndex = -1
    LineToCellIntersection.SurfaceIndex = -1
    For SurfaceIndex = 0 To C.NumberOfSurfaces - 1
        Dim ASurfacePlane As Plane
        ASurfacePlane = ApproximateFractionPlane(C.Surfaces(SurfaceIndex).SurfacePlane)
        If PointToSurfaceDistance(L.StartPoint, ASurfacePlane) > 0 Then 'Back face culling
            Hit = LineToPlaneIntersection(L, ASurfacePlane) 'Ray tracing
            If Hit.Collided Then
                Dim Inside As Boolean
                Inside = True
                For S2 = 0 To C.NumberOfSurfaces - 1 'Clipping
                    If S2 <> SurfaceIndex And PointToSurfaceDistance(Hit.PointOfImpact, ApproximateFractionPlane(C.Surfaces(S2).SurfacePlane)) > 0 Then
                        Inside = False
                        Exit For
                    End If
                Next S2
                If Inside Then
                    LineToCellIntersection.CellIndex = CellIndex
                    LineToCellIntersection.SurfaceIndex = SurfaceIndex
                    LineToCellIntersection.PlaneIntersection = Hit
                    L.EndPoint = Hit.PointOfImpact
                End If
            End If
        End If
    Next SurfaceIndex
End Function

Public Function RemovePlaneFromCell(C As Cell, SurfaceIndex As Long) As Cell
    RemovePlaneFromCell = C
    If SurfaceIndex < RemovePlaneFromCell.NumberOfSurfaces - 1 Then
        RemovePlaneFromCell.Surfaces(SurfaceIndex) = RemovePlaneFromCell.Surfaces(RemovePlaneFromCell.NumberOfSurfaces - 1)
    End If
    RemovePlaneFromCell.NumberOfSurfaces = RemovePlaneFromCell.NumberOfSurfaces - 1
    RemovePlaneFromCell = GenerateCellFromPlanes(RemovePlaneFromCell)
    UpdateCellBound RemovePlaneFromCell
End Function

'If the edges were made of a metal wire around the polyogn then the result would be its center of mass
Public Function GetSurfaceCenter(C As Cell, S As Surface) As Vector3
    Dim P As Long
    Dim P2 As Long
    Dim TotalWeight As Double
    Dim Sum As Vector3
    If S.NumberOfIndices = 0 Then
        GetSurfaceCenter = MakeVector3(0, 0, 0)
        Exit Function
    End If
    For P = 0 To S.NumberOfIndices - 1
        Dim StartPos As Vector3
        Dim EndPos As Vector3
        Dim Length As Double
        StartPos = C.Positions(S.Indices(P))
        EndPos = C.Positions(S.Indices((P + 1) Mod S.NumberOfIndices))
        Length = AbsVector3(SubVector3(EndPos, StartPos))
        TotalWeight = TotalWeight + Length
        Sum = AddVector3(Sum, MulVector3(MiddleVector3(StartPos, EndPos), Length))
    Next P
    GetSurfaceCenter = DivVector3(Sum, TotalWeight)
End Function

Public Function AddPointToSurface(ByRef C As Cell, ByRef S As Surface, Position As Vector3) As Long
    AddPointToSurface = AddPosition(C, Position)
    AddPointToSurfaceByIndex S, AddPointToSurface
End Function

Public Sub AddPointToSurfaceByIndex(ByRef S As Surface, Index As Long)
    Dim I As Long
    'Remove duplicates in case that crossing a point adds the same index multiple times from different edge-plane intersections
    For I = 0 To S.NumberOfIndices - 1
        If S.Indices(I) = Index Then
            Exit Sub
        End If
    Next I
    S.Indices(S.NumberOfIndices) = Index
    S.NumberOfIndices = S.NumberOfIndices + 1
End Sub

'Precondition: CellIntersectsPlane(SourceCell, CuttingPlane) = true to avoid redundant planes
'Postcondition: Returns the volume from SourceCell that is on the positive side of CuttingPlane but not from the negative side
Public Function GetFragmentOfCell(SourceCell As Cell, CuttingPlane As FractionPlane) As Cell
    Dim SurfaceIndex As Long
    Dim WS As Long
    Dim P As Long
    Dim Dist As Double
    Dim NextDist As Double
    Dim Pos As Vector3
    Dim NextPos As Vector3
    Dim CutPos As Vector3
    Dim HasFront As Boolean
    Dim HasBack As Boolean
    Dim Hole As HoleFillerData
    Dim ACuttingPlane As Plane
    ACuttingPlane = ApproximateFractionPlane(CuttingPlane)
    GetFragmentOfCell.Settings = SourceCell.Settings
    
    Debug.Assert CellIntersectsPlane(SourceCell, CuttingPlane)
    
    For SurfaceIndex = 0 To SourceCell.NumberOfSurfaces - 1
        HasFront = False
        HasBack = False
        For P = 0 To SourceCell.Surfaces(SurfaceIndex).NumberOfIndices - 1
            Pos = SourceCell.Positions(SourceCell.Surfaces(SurfaceIndex).Indices(P))
            Dist = PointToSurfaceDistance(Pos, ACuttingPlane)
            If Dist > Geometry_OnPlaneTolerance Then
                HasFront = True
            ElseIf Dist < -Geometry_OnPlaneTolerance Then
                HasBack = True
            End If
        Next P
        If HasFront And HasBack Then 'Cut the surface along the plane so that the part on the negative side of the plane is removed
            WS = GetFragmentOfCell.NumberOfSurfaces
            GetFragmentOfCell.NumberOfSurfaces = GetFragmentOfCell.NumberOfSurfaces + 1
            'Copy the surface plane
            GetFragmentOfCell.Surfaces(WS).SurfacePlane = SourceCell.Surfaces(SurfaceIndex).SurfacePlane
            GetFragmentOfCell.Surfaces(WS).SurfaceType = SourceCell.Surfaces(SurfaceIndex).SurfaceType
            'Copy the vertices
            For P = 0 To SourceCell.Surfaces(SurfaceIndex).NumberOfIndices - 1
                Pos = SourceCell.Positions(SourceCell.Surfaces(SurfaceIndex).Indices(P))
                NextPos = SourceCell.Positions(SourceCell.Surfaces(SurfaceIndex).Indices((P + 1) Mod SourceCell.Surfaces(SurfaceIndex).NumberOfIndices))
                Dist = PointToSurfaceDistance(Pos, ACuttingPlane)
                NextDist = PointToSurfaceDistance(NextPos, ACuttingPlane)
                If Dist < 0 And NextDist >= 0 Then
                    'Cut along an edge
                    CutPos = LerpVector3(Pos, NextPos, Dist / (Dist - NextDist))
                    HoleFiller_AddCut Hole, AddPointToSurface(GetFragmentOfCell, GetFragmentOfCell.Surfaces(WS), CutPos)
                End If
                If Dist >= 0 Then
                    'Whole vertex
                    AddPointToSurface GetFragmentOfCell, GetFragmentOfCell.Surfaces(WS), Pos
                    If NextDist < 0 Then
                        'Cut along an edge
                        CutPos = LerpVector3(Pos, NextPos, Dist / (Dist - NextDist))
                        HoleFiller_AddCut Hole, AddPointToSurface(GetFragmentOfCell, GetFragmentOfCell.Surfaces(WS), CutPos)
                    End If
                End If
            Next P
        ElseIf Not HasBack Then 'Copy the whole surface while assigning new indices
            WS = GetFragmentOfCell.NumberOfSurfaces
            GetFragmentOfCell.NumberOfSurfaces = GetFragmentOfCell.NumberOfSurfaces + 1
            'Copy the surface
            GetFragmentOfCell.Surfaces(WS).SurfacePlane = SourceCell.Surfaces(SurfaceIndex).SurfacePlane
            GetFragmentOfCell.Surfaces(WS).SurfaceType = SourceCell.Surfaces(SurfaceIndex).SurfaceType
            'Copy the vertices
            For P = 0 To SourceCell.Surfaces(SurfaceIndex).NumberOfIndices - 1
                'Add whole vertex
                Pos = SourceCell.Positions(SourceCell.Surfaces(SurfaceIndex).Indices(P))
                Dim TargetPointIndex As Long
                TargetPointIndex = AddPointToSurface(GetFragmentOfCell, GetFragmentOfCell.Surfaces(WS), Pos)
                'Find edges defining the hole
                Dist = PointToSurfaceDistance(Pos, ACuttingPlane)
                If Abs(Dist) < Geometry_WeldTolerance Then
                    'Cut into a corner
                    Dim Seed As Integer
                    Seed = Int(Rnd * 10000)
                    HoleFiller_AddCut Hole, TargetPointIndex
                End If
            Next P
        End If
    Next SurfaceIndex
    'Use the cut indices collected in Hole to fill the hole
    HoleFiller_Finalize Hole, GetFragmentOfCell, CuttingPlane
    Dim Message As String
    Message = ""
    If Not ValidateCell(GetFragmentOfCell, Message) Then
        MsgBox Message, vbCritical, "Validation of cell fragment failed!"
        GetFragmentOfCell.Settings.IsBroken = True
    End If
    UpdateCellBound GetFragmentOfCell
End Function


'Returns true iff A is inside of B
Public Function CellInsideCell(A As Cell, B As Cell) As Boolean
    CellInsideCell = True
    Dim P As Long
    For P = 0 To A.NumberOfPositions - 1
        If Not PointInsideCell(A.Positions(P), B) Then
            CellInsideCell = False
            Exit Function
        End If
    Next P
End Function

'Returns the combintation of CellA And CellB
Public Function CombineCells(CellA As Cell, CellB As Cell) As Cell
    Dim CurrentCell As Cell
    Dim ASurface As Plane
    CurrentCell = MakeCell_Box(MakeFractionBox_Big, CellA.Settings, True)
    Dim SurfaceIndex As Long
    For SurfaceIndex = 0 To CellA.NumberOfSurfaces - 1
        If CellA.Surfaces(SurfaceIndex).SurfaceType <> Geometry_SurfaceType_Infinity And CellOnPlaneNegative(CellB, CellA.Surfaces(SurfaceIndex).SurfacePlane) And CellIntersectsPlane(CurrentCell, CellA.Surfaces(SurfaceIndex).SurfacePlane) Then
            CurrentCell = GetFragmentOfCell(CurrentCell, FlipFractionPlane(CellA.Surfaces(SurfaceIndex).SurfacePlane))
        End If
    Next SurfaceIndex
    For SurfaceIndex = 0 To CellA.NumberOfSurfaces - 1
        If CellB.Surfaces(SurfaceIndex).SurfaceType <> Geometry_SurfaceType_Infinity And CellOnPlaneNegative(CellA, CellB.Surfaces(SurfaceIndex).SurfacePlane) And CellIntersectsPlane(CurrentCell, CellB.Surfaces(SurfaceIndex).SurfacePlane) Then
            CurrentCell = GetFragmentOfCell(CurrentCell, FlipFractionPlane(CellB.Surfaces(SurfaceIndex).SurfacePlane))
        End If
    Next SurfaceIndex
    CombineCells = CurrentCell
    UpdateCellBound CombineCells
    'いい Make a sanity check on the bounding boxes so that going to infinity will raise an exception
    '   The cells must have a polygon in common and have no overlap in volume
End Function

'Returns the only plane in CellA that matches a negation of a plane in CellB
'Returns the index of the only surface in CellA with the whole CellB on the positive side or -1 if none or multiple
Public Function FindGapIndexA(CellA As Cell, CellB As Cell) As Long
    Dim SurfaceIndex As Long
    FindGapIndexA = -1 'No result
    For SurfaceIndex = 0 To CellA.NumberOfSurfaces - 1
        If CellTouchesPlane(CellB, CellA.Surfaces(SurfaceIndex).SurfacePlane) Then
            If FindGapIndexA > -1 Then
                FindGapIndexA = -1 'Multiple results
                Exit Function
            Else
                FindGapIndexA = SurfaceIndex 'First result
            End If
        End If
    Next SurfaceIndex
End Function

Public Function GetTranslationOfCell(C As Cell, Translation As FractionPoint) As Cell
    GetTranslationOfCell = C
    Dim S As Long
    'Translate each plane
    For S = 0 To C.NumberOfSurfaces - 1
        GetTranslationOfCell.Surfaces(S).SurfacePlane.Normal = C.Surfaces(S).SurfacePlane.Normal
        GetTranslationOfCell.Surfaces(S).SurfacePlane.Point = AddFractionPoint(C.Surfaces(S).SurfacePlane.Point, Translation)
    Next S
    GetTranslationOfCell = GenerateCellFromPlanes(GetTranslationOfCell)
    UpdateCellBound GetTranslationOfCell
    Dim Message As String
    Message = ""
    If Not ValidateCell(GetTranslationOfCell, Message) Then
        MsgBox Message, vbCritical, "Validation of cell transform failed!"
        GetTranslationOfCell.Settings.IsBroken = True
    End If
End Function

Public Function GetReflectionOfCell(C As Cell, Center As FractionPoint, Dimension As Integer) As Cell
    GetReflectionOfCell = C
    Dim SurfaceIndex As Long
    For SurfaceIndex = 0 To C.NumberOfSurfaces - 1
        GetReflectionOfCell.Surfaces(SurfaceIndex).SurfacePlane = MirrorFractionPlane(C.Surfaces(SurfaceIndex).SurfacePlane, Center, Dimension)
    Next SurfaceIndex
    GetReflectionOfCell = GenerateCellFromPlanes(GetReflectionOfCell)
    UpdateCellBound GetReflectionOfCell
    Dim Message As String
    Message = ""
    If Not ValidateCell(GetReflectionOfCell, Message) Then
        MsgBox Message, vbCritical, "Validation of cell transform failed!"
        GetReflectionOfCell.Settings.IsBroken = True
    End If
End Function
