単独な3D曲線の取得する3

明けましておめでとうございます。 ﾏｲﾍﾟｰｽでやっていきます。
こちらの続きです。
単独な3D曲線の取得する2 - C#ATIA

昨年末に線形8分木を実装しましたが、若干効率が悪いことには
気が付いていた為、処理効率を上げるように修正してみました。

CATMainを持つ標準ﾓｼﾞｭｰﾙ "FindSingle3DCurve.bas" です。

'vba FindSingle3DCurve Ver0.0.3  using-'KCL0.09' 'FindSingle3DCurve_Octree0.0.3'
'単独な3D曲線の取得

Option Explicit

'*** 設定 ***
Private Const Tolerance = 0.001                         '端点一致ﾄﾚﾗﾝｽ
Private Const InSpaceMaxCount = 200&                    '同一空間内最大数(目安)
Private Const SingleCurveFooter = " - SingleCurve!!"    '単独曲線名に追記する文字
'************

Sub CATMain()
    'ﾄﾞｷｭﾒﾝﾄのﾁｪｯｸ
    If Not CanExecute(Array("PartDocument", "ProductDocument")) Then Exit Sub
    Dim Msg$
    
    '形状ｾｯﾄの選択
    Msg = "3D曲線の入っている形状ｾｯﾄを指定してください : ESCｷｰ 終了"
    Dim Hb As HybridBody
    Set Hb = KCL.SelectItem(Msg, "HybridBody")
    If KCL.IsNothing(Hb) Then Exit Sub
    
    '線の取得
    Dim Lines As Collection: Set Lines = GetLines(Hb)
    If KCL.IsNothing(Lines) Then
        Msg = "指定した' " + Hb.Name + " 'には、3D曲線がありませんでした!"
        MsgBox Msg: Exit Sub
    End If
    
    '警告
    Msg = CStr(Lines.Count) + "本あります。実行しますか?"
    
    If MsgBox(Msg, vbYesNo + vbInformation) = vbNo Then Exit Sub
    Debug.Print "** start003 ** : " + CStr(Lines.Count) + "本" 'debug
    Debug.Print "InSpaceMaxCount : " + CStr(InSpaceMaxCount) 'debug
    KCL.SW_Start 'debug
    
    '端点座標取得
    Dim EndPnts As Collection: Set EndPnts = GetEndPntCoordinates(Lines)
    Debug.Print "EndPnts : " + CStr(KCL.SW_GetTime) 'debug
    
    '組み合わせﾘｽﾄ取得
    Dim OctList As Collection
    Set OctList = FindSingle3DCurve_Octree.GetLinerOctreeList(EndPnts, Tolerance, InSpaceMaxCount)
    If KCL.IsNothing(OctList) Then
        Msg = "ｴﾗｰです。ﾏｸﾛ内の設定値を修正して下さい(特に 'OctLv' )"
        MsgBox Msg: Exit Sub
    End If
    Debug.Print "OctList : " + CStr(KCL.SW_GetTime) 'debug
    
    '単独曲線ｲﾝﾃﾞｯｸｽ取得
    Dim SingleList As Collection: Set SingleList = GetSingleIdx(EndPnts, OctList)
    Debug.Print "SingleList : " + CStr(SingleList.Count) + " : " + CStr(KCL.SW_GetTime)
    If SingleList.Count < 1 Then
        Msg = "指定した' " + Hb.Name + " 'には、単独な3D曲線が見つかりませんでした!" + vbNewLine + _
              "(ﾄﾚﾗﾝｽ " & Tolerance & " mm)"
        MsgBox Msg: Exit Sub
    End If
    
    '確認
    Msg = SingleList.Count & "本の単独な3D曲線が見つかりました!!" + vbNewLine + _
          "名前の最後に ' " + SingleCurveFooter + " ' を追記しますか?"
    If MsgBox(Msg, vbYesNo + vbInformation) = vbNo Then GoTo MacroEnd
    
    'ﾘﾈｰﾑ
    Call AppendName(Lines, SingleList)
    
MacroEnd:
    'Call Dump(ByVal SingleList, ByVal Lines, ByVal EndPnts) 'debug
    'Call DumpCount(ByVal OctList) 'debug
    MsgBox "終了!"
End Sub

'単独曲線のﾘﾈｰﾑ
Private Sub AppendName(ByVal Lines As Collection, ByVal IdxList As Collection)
    Dim i
    For Each i In IdxList
        With Lines.Item(i)
            .Name = .Name & SingleCurveFooter
        End With
    Next
End Sub

'初期化済み配列生成 - ｵﾌﾞｼﾞｪｸﾄNG
Private Function InitRangeAry(ByVal Count&, ByVal Value As Variant)
    Dim ary() As Variant: ReDim ary(Count)
    Dim i&
    For i = 0 To Count
        ary(i) = Value
    Next
    InitRangeAry = ary
End Function

'単独曲線端点とそうではない端点を分ける
''' @return:Collection(Collection-一致,Collection-未確定)
Private Function GroupBySinglePointList(ByVal PntList As Collection, _
                                        ByVal EndPnts As Collection, _
                                        ByVal CombAry) As Variant
    Dim TrueList As Collection: Set TrueList = New Collection
    Dim FalseList As Collection: Set FalseList = New Collection
    Dim Idx
    For Each Idx In PntList
        If CombAry(EndPnts.Item(Idx)(3)) Then
            TrueList.Add Idx
        Else
            FalseList.Add Idx
        End If
    Next
    GroupBySinglePointList = Array(TrueList, FalseList)
End Function

'単独曲線ｲﾝﾃﾞｯｸｽ取得
Private Function GetSingleIdx(ByVal EndPnts As Collection, ByVal OctList As Collection) As Collection
    Dim i&, j&, i_crv_idx, j_crv_idx
    Dim CombAry As Variant: CombAry = InitRangeAry(EndPnts.Count * 0.5, False)
    Dim Space As Collection
    Dim GroupAry, TrueList As Collection, FalseList As Collection
    
    For Each Space In OctList
        '未確定同士
        GroupAry = GroupBySinglePointList(Space, EndPnts, CombAry)
        Set TrueList = GroupAry(0)
        Set FalseList = GroupAry(1)
        For i = 1 To FalseList.Count
            i_crv_idx = EndPnts.Item(FalseList.Item(i))(3)
            For j = i + 1 To FalseList.Count
                j_crv_idx = EndPnts.Item(FalseList.Item(j))(3)
                If i_crv_idx = j_crv_idx Then GoTo Continue1
                If CombAry(i_crv_idx) And CombAry(j_crv_idx) Then GoTo Continue1
                If IsPosEqual(EndPnts.Item(FalseList.Item(i)), EndPnts.Item(FalseList.Item(j))) Then
                    CombAry(i_crv_idx) = True
                    CombAry(j_crv_idx) = True
                End If
Continue1:
            Next
        Next
        
        '未確定-確定
        GroupAry = GroupBySinglePointList(FalseList, EndPnts, CombAry)
        Set FalseList = GroupAry(1)
        For i = 1 To FalseList.Count
            i_crv_idx = EndPnts.Item(FalseList.Item(i))(3)
            For j = 1 To TrueList.Count
                j_crv_idx = EndPnts.Item(TrueList.Item(j))(3)
                If i_crv_idx = j_crv_idx Then GoTo Continue2
                If CombAry(i_crv_idx) And CombAry(j_crv_idx) Then GoTo Continue2
                If IsPosEqual(EndPnts.Item(FalseList.Item(i)), EndPnts.Item(TrueList.Item(j))) Then
                    CombAry(i_crv_idx) = True
                End If
Continue2:
            Next
        Next
    Next
    
    Dim SingleIdxList As Collection: Set SingleIdxList = New Collection
    For i = 1 To UBound(CombAry)
        If CombAry(i) = False Then SingleIdxList.Add i
    Next
    Set GetSingleIdx = SingleIdxList
End Function

'座標値一致
Private Function IsPosEqual(ByVal p1 As Variant, ByVal p2 As Variant) As Boolean
    IsPosEqual = False
    Dim i&
    For i = 0 To 2
        If Abs(p2(i) - p1(i)) > Tolerance Then Exit Function
    Next
    If Sqr(Abs((p2(0) - p1(0)) * (p2(0) - p1(0)) + _
               (p2(1) - p1(1)) * (p2(1) - p1(1)) + _
               (p2(2) - p1(2)) * (p2(2) - p1(2)))) < Tolerance Then
        IsPosEqual = True
    End If
End Function

'端点取得
Private Function GetEndPntCoordinates(ByRef Lines As Collection) As Collection
    Dim Doc As PartDocument: Set Doc = KCL.GetParent_Of_T(Lines.Item(1), "PartDocument")
    Dim SPAWB As Workbench: Set SPAWB = Doc.GetWorkbench("SPAWorkbench")
    Dim EndPnt As Collection: Set EndPnt = New Collection
    Dim Pos(8) As Variant
    Dim i&
    For i = 1 To Lines.Count
        Call SPAWB.GetMeasurable(Lines.Item(i)).GetPointsOnCurve(Pos)
        EndPnt.Add KCL.JoinAry(KCL.GetRangeAry(Pos, 0, 2), Array(i))
        EndPnt.Add KCL.JoinAry(KCL.GetRangeAry(Pos, 6, 8), Array(i))
    Next
    Set GetEndPntCoordinates = EndPnt
End Function

'線取得
Private Function GetLines(ByVal Hb As HybridBody) As Collection
    Dim Hss As HybridShapes: Set Hss = Hb.HybridShapes
    If Hss.Count < 1 Then Exit Function
    
    Dim Pt As Part: Set Pt = KCL.GetParent_Of_T(Hb, "PartDocument").Part
    Dim Fact As HybridShapeFactory: Set Fact = Pt.HybridShapeFactory
    Dim Lines As Collection: Set Lines = New Collection
    Dim Hs As HybridShape
    For Each Hs In Hss
        Select Case Fact.GetGeometricalFeatureType(Hs)
            Case 2, 3, 4
                Lines.Add Hs
        End Select
    Next
    Set GetLines = Lines
End Function

'debug用
Private Sub Dump(ByVal SingleList, ByVal Lines, ByVal EndPnts)
    Dim ss
    For Each ss In SingleList
        Debug.Print ss & ":" & Lines.Item(EndPnts.Item(ss)(3)).Name
    Next
End Sub

'debug用
Private Sub DumpCount(ByVal List)
    Dim ss
    For Each ss In List
        Debug.Print ss.Count
    Next
End Sub

'ver0.0.1 - 総当りのﾀﾞﾒなやつ
'ver0.0.2 - 線形8分木を実装
'ver0.0.3 - 端点一致ｱﾙｺﾞﾘｽﾞﾑ効率化

続いて、組み合わせﾘｽﾄを作り出す "FindSingle3DCurve_Octree.bas" です。

'vba FindSingle3DCurve_Octree Ver0.0.3  using-'KCL0.09'
'ﾓｰﾄﾝ順序を利用した8分木空間分割ﾗｲﾌﾞﾗﾘ
'FindSingle3DCurveﾏｸﾛ専用です

Option Explicit

Private Const CLINER8TREEMANAGER_MAXLEVEL = 7   '有効空間分割最大ﾚﾍﾞﾙ

Private m_Level&                                '分割ﾚﾍﾞﾙ
Private m_Tolerance#                            '一致ﾄﾚﾗﾝｽ
Private m_MaxCount&                             '同一空間内最大数(目安)
Private m_AxisCount&                            '空間分割時の各軸の最大数
Private m_MinPos                                '空間最小座標
Private m_Unit                                  '空間単位ｻｲｽﾞ
Private m_ToleranceRatio                        '空間単位ｻｲｽﾞに対してのﾄﾚﾗﾝｽ比率
Private m_CellCount&()                          'ﾚﾍﾞﾙ毎の空間数

'*** Octree ***
'http://marupeke296.com/COL_3D_No15_Octree.html

'線形8分木ﾘｽﾄ取得
''' @param :Pnts-Collection(array(Double))-座標値郡
''' @param :Tolerance-Double-一致ﾄﾚﾗﾝｽ
''' @param :MaxCount-long-同一空間内最大数(目安)
''' @return:Collection(Collection(long))-空間別点Idx郡
Function GetLinerOctreeList(ByVal Pnts As Collection, ByVal Tolerance#, ByVal MaxCount&)
    Set GetLinerOctreeList = Nothing
    
    '座標値郡のIdxList作成
    Dim PntIdxList As Collection: Set PntIdxList = InitRangeList(Pnts.Count)
    
    '空間内最大数以下の場合、そのまま返す
    Dim DecidedList As Collection: Set DecidedList = New Collection
    If Pnts.Count < MaxCount Then
        Call DecidedList.Add(PntIdxList)
        GoTo FuncEnd
    End If
    
    '初期設定
    If Not SetStart(Tolerance, MaxCount) Then
        MsgBox "設定値が不正です"
        Exit Function
    End If
    
    'IdxListをﾙｰﾄ空間として登録
    Dim CheckList As Collection: Set CheckList = New Collection
    Call CheckList.Add(PntIdxList)
    
    'ﾚﾍﾞﾙの応じたIdx用配列
    Dim SpeceEnum, TempList As Collection, ReCheckList As Collection
    Dim Space, TempSpace, SpaAry, i&, Idx
    Do
        Set ReCheckList = New Collection
        '配置
        For Each Space In CheckList
            If Not SetSpaceInfo(Pnts, Space) Then '空間が小さすぎる
                Call DecidedList.Add(Space)
                GoTo Continue
            End If
            
            SpeceEnum = InitRangeAry(m_CellCount(m_Level), -1)
            Set TempList = New Collection
            
            For Each Idx In Space
                SpaAry = GetMortonNum(Pnts(Idx))
                For i = 0 To UBound(SpaAry)
                    If SpeceEnum(SpaAry(i)) < 0 Then
                        Call TempList.Add(InitSpace())
                        SpeceEnum(SpaAry(i)) = TempList.Count
                    End If
                    Call TempList.Item(SpeceEnum(SpaAry(i))).Add(Idx)
                Next
            Next
            
            '空間毎の数が多いものを再配置
            For Each TempSpace In TempList
                Select Case True
                    Case TempSpace.Count < 2
                        '空間に1個しかないものは検証しない
                        'Call DecidedList.add(TempSpace)
                    Case TempSpace.Count > m_MaxCount
                        Call ReCheckList.Add(TempSpace)
                    Case Else
                        Call DecidedList.Add(TempSpace)
                End Select
            Next
Continue:
        Next
        If ReCheckList.Count < 1 Then Exit Do
        Set CheckList = ReCheckList
    Loop
    
FuncEnd:
    Call Q_ISort_List(DecidedList)
    Set GetLinerOctreeList = DecidedList
    Set PntIdxList = Nothing
    Set DecidedList = Nothing
    Set CheckList = Nothing
    Set ReCheckList = Nothing
    Set TempList = Nothing
    ReDim SpeceEnum(0)
End Function

'空のｺﾚｸｼｮﾝ
Private Function InitSpace() As Collection
    Set InitSpace = New Collection
End Function

'線形8分木準備
''' @param :Level-long-分割ﾚﾍﾞﾙ
''' @param :Tolerance-Double-一致ﾄﾚﾗﾝｽ
''' @param :MaxCount-long-同一空間内最大数(目安)
''' @return:Boolean
Private Function SetStart(ByVal Tolerance#, ByVal MaxCount&) As Boolean
    SetStart = False
    If Tolerance <= 0 Then Exit Function
    
    m_Tolerance = Tolerance
    m_MaxCount = MaxCount
    m_ToleranceRatio = InitRangeAry(2, 0)
    
    ReDim m_CellCount(CLINER8TREEMANAGER_MAXLEVEL + 1)
    m_CellCount(0) = 1
    Dim i&
    For i = 1 To UBound(m_CellCount)
        m_CellCount(i) = m_CellCount(i - 1) * 8
    Next
    SetStart = True
End Function

'空間ｻｲｽﾞとﾄﾚﾗﾝｽからﾚﾍﾞﾙ算出し設定
Private Function SetLevel(ByVal W) As Boolean
    SetLevel = False
    Dim Min#: Min = 1.79769313486231E+308
    Dim i&
    For i = 0 To 2
         If Min > W(i) Then Min = W(i)
    Next
    Dim TmpLv&: TmpLv = Fix(Log_n((Min / (m_Tolerance * 2 + 0.002)), 2))
    If TmpLv > CLINER8TREEMANAGER_MAXLEVEL Then
        m_Level = CLINER8TREEMANAGER_MAXLEVEL
    Else
        m_Level = TmpLv
    End If
    
    If m_Level < 1 Then Exit Function
    m_AxisCount = sl(1, m_Level)
    SetLevel = True
End Function

'空間情報設定
''' @param :Pnts-Collection(array(Double))-座標値郡
''' @param :Idxs-Collection(long)-座標値郡Idx
''' @return:Boolean
Private Function SetSpaceInfo(ByVal Pnts, ByVal Idxs As Collection) As Boolean
    SetSpaceInfo = False
    Dim SpSize: SpSize = GetSpaceSize_Idx(Pnts, Idxs)
    m_MinPos = AryAdd(SpSize(0), m_Tolerance * -1)
    Dim W: W = ArySub(AryAdd(SpSize(1), m_Tolerance), m_MinPos)
    If Not SetLevel(W) Then Exit Function
    m_Unit = AryDiv(W, m_AxisCount)
    
    Dim i&
    For i = 0 To 2
        m_ToleranceRatio(i) = m_ToleranceRatio(i) / m_Unit(i)
    Next
    SetSpaceInfo = True
End Function

'座標→最小ﾚﾍﾞﾙ空間番号郡取得
''' @param :Pos-array(Double)-座標値
''' @return:array(Long)
Private Function GetMortonNum(ByVal Pos As Variant) As Variant
    Dim Ratio#(2), Inte&(2), Dec#(2)
    Dim i&
    For i = 0 To 2
        Ratio(i) = (Pos(i) - m_MinPos(i)) / m_Unit(i)
        Inte(i) = Fix(Ratio(i))
        Dec(i) = Ratio(i) - Inte(i)
    Next
    
    Dim Axis(2) As Variant
    Dim AxisNums As Collection
    For i = 0 To 2
        Set AxisNums = New Collection
        AxisNums.Add Inte(i)
        If Dec(i) <= m_ToleranceRatio(i) And Inte(i) > 0 Then AxisNums.Add Inte(i) - 1
        If 1# - Dec(i) <= m_ToleranceRatio(i) And Inte(i) < m_AxisCount - 1 Then AxisNums.Add Inte(i) + 1
        Set Axis(i) = AxisNums
    Next
    
    Dim x, y, z
    Dim SpaNo(): ReDim SpaNo(Axis(0).Count * Axis(1).Count * Axis(2).Count)
    Dim cnt&: cnt = -1
    For Each x In Axis(0)
        For Each y In Axis(1)
            For Each z In Axis(2)
                cnt = cnt + 1
                SpaNo(cnt) = Get3DMortonNumber(x, y, z)
                If SpaNo(cnt) < 0 And SpaNo(cnt) >= m_CellCount(m_Level) Then cnt = cnt - 1
    Next: Next: Next
    ReDim Preserve SpaNo(cnt)
    GetMortonNum = SpaNo
End Function

'ビット分割関数
''' @param :n-long
''' @return:long
Private Function BitSeparateFor3D(ByVal n&) As Long
    Dim s As Long: s = n
    s = (s Or sl(s, 8)) And &HF00F
    s = (s Or sl(s, 4)) And &HC30C3
    s = (s Or sl(s, 2)) And &H249249
    BitSeparateFor3D = s
End Function

'8分木モートン順序算出関数
''' @param :x-long
''' @param :y-long
''' @param :z-long
''' @return:long
Private Function Get3DMortonNumber(ByVal x&, ByVal y&, ByVal z&) As Long
   Get3DMortonNumber = BitSeparateFor3D(x) Or _
                       sl(BitSeparateFor3D(y), 1) Or _
                       sl(BitSeparateFor3D(z), 2)
End Function

'座標値郡から空間ｻｲｽﾞ取得
''' @param :EndPnts-Collection(array(Double))-座標値郡
''' @return:array(array(Double))-0:最小値 1:最大値
Private Function GetSpaceSize_Idx(ByVal EndPnts As Collection, ByVal IdxList As Collection) As Variant
    Dim Min: Min = InitRangeAry(2, 1.79769313486231E+308)
    Dim Max: Max = InitRangeAry(2, -1.79769313486231E+308)
    
    Dim Idx, i&
    For Each Idx In IdxList
        For i = 0 To 2
            If Min(i) > EndPnts.Item(Idx)(i) Then Min(i) = EndPnts.Item(Idx)(i)
            If Max(i) < EndPnts.Item(Idx)(i) Then Max(i) = EndPnts.Item(Idx)(i)
        Next
    Next
    GetSpaceSize_Idx = Array(Min, Max)
End Function

'*** BitShift ***
'http://www.geocities.co.jp/SiliconValley/4334/unibon/asp/bitshift2.html
' 左シフト
Private Function sl(ByVal x&, ByVal n&) As Long
    If n = 0 Then
        sl = x
    Else
        Dim k: k = CLng(2 ^ (32 - n - 1))
        Dim d: d = x And (k - 1)
        Dim c: c = d * CLng(2 ^ n)
        If x And k Then c = c Or &H80000000
        sl = c
    End If
End Function

'*** VBA不足関数 ***
'配列同士の引き算-細かいチェック無し
Private Function ArySub(ByVal A, ByVal B) As Variant
    ArySub = Array(A(0) - B(0), A(1) - B(1), A(2) - B(2))
End Function

'配列と実数の足し算-細かいチェック無し
Private Function AryAdd(ByVal A, ByVal B#) As Variant
    AryAdd = Array(A(0) + B, A(1) + B, A(2) + B)
End Function

'配列と実数の割り算-細かいチェック無し
Private Function AryDiv(ByVal A, ByVal B#) As Variant
    AryDiv = Array(A(0) / B, A(1) / B, A(2) / B)
End Function

'初期化済み配列生成 - ｵﾌﾞｼﾞｪｸﾄNG
Private Function InitRangeAry(ByVal Count&, ByVal Value As Variant)
    Dim ary() As Variant: ReDim ary(Count)
    Dim i&
        For i = 0 To Count
            ary(i) = Value
        Next
    InitRangeAry = ary
End Function

'初期化済みｺﾚｸｼｮﾝ生成
Private Function InitRangeList(ByVal Count&) As Collection
    Dim List As Collection: Set List = New Collection
    Dim i&
    For i = 1 To Count
        List.Add i
    Next
    Set InitRangeList = List
End Function

'ｺﾚｸｼｮﾝｶｳﾝﾀによるｸｲｯｸｿｰﾄ　非再帰版
Private Sub Q_ISort_List(ByRef List As Collection)
    Dim Threashold&: Threashold = 64
    Dim Stack As Collection: Set Stack = New Collection
    Stack.Add 1, CStr(Stack.Count + 1)
    Stack.Add List.Count, CStr(Stack.Count + 1)
    
    Dim Pivot, Temp1, Temp2
    Dim LeftIdx&, RightIdx&, i&, j&
    Do While Stack.Count > 0
        LeftIdx = Stack(CStr(Stack.Count - 1))
        RightIdx = Stack(CStr(Stack.Count))
        Stack.Remove Stack.Count
        Stack.Remove Stack.Count
        'クイックソート
        If LeftIdx < RightIdx Then
            Set Pivot = List((LeftIdx + RightIdx) / 2)
            i = LeftIdx
            j = RightIdx
            
            Do While i <= j
                Do While List(i).Count < Pivot.Count
                    i = i + 1
                Loop
                Do While List(j).Count > Pivot.Count
                    j = j - 1
                Loop
                If i <= j Then
                    Set Temp1 = List(i)
                    Set Temp2 = List(j)
                    List.Add Temp1, After:=j
                    List.Remove j
                    List.Add Temp2, After:=i
                    List.Remove i
                    i = i + 1
                    j = j - 1
                End If
            Loop
            
            If RightIdx - i >= 0 Then
                If RightIdx - i <= Threashold Then
                    ComboInsertionSort List, i, RightIdx
                Else
                    Stack.Add i, CStr(Stack.Count + 1)
                    Stack.Add RightIdx, CStr(Stack.Count + 1)
                End If
            End If
            
            If j - LeftIdx >= 0 Then
                If j - LeftIdx <= Threashold Then
                    ComboInsertionSort List, LeftIdx, j
                Else
                    Stack.Add LeftIdx, CStr(Stack.Count + 1)
                    Stack.Add j, CStr(Stack.Count + 1)
                End If
            End If
        End If
    Loop
End Sub

'InsertSort
Private Sub ComboInsertionSort(ByRef List, ByVal MinIdx&, ByVal MaxIdx&)
    Dim Temp1, Temp2
    Dim i&, j&: j = 1
    For j = MinIdx To MaxIdx
        i = j - 1
        Do While i >= 1
            If List(i + 1).Count < List(i).Count Then
                Set Temp1 = List(i + 1)
                Set Temp2 = List(i)
                List.Add Temp2, After:=i + 1
                List.Remove i + 1
                List.Add Temp1, After:=i
                List.Remove i
            Else
                Exit Do
            End If
            i = i - 1
        Loop
    Next
End Sub

'nを底とする対数
Private Function Log_n(x, n)
    Log_n = Log(x) / Log(n)
End Function

'ver0.0.1 - とりあえず完成
'ver0.0.3 - 分割ﾚﾍﾞﾙ自動化,線形8分木ｿｰﾄ化

細々した事を書いても仕方ない為説明を割愛しますが、Ver0.0.2の際は
線形8分木ﾘｽﾄを取得する際、事前にInit関数を利用し分割ﾚﾍﾞﾙを指定する
必要がありましたが、与えられたﾃﾞｰﾀ(端点座標値郡)と一致ﾄﾚﾗﾝｽから
自動的に可能な限り深いﾚﾍﾞﾙまで分割するようにしました。
(そもそもｲﾝｽﾀﾝｽを作成しないのに、" Init " と言う名称の関数がある事
自体おかしかったです…)

Ver0.0.2と0.0.3を比べた結果はこちらです。

Ver0.0.2
** start ** : 4100本
OctLv : 7  InSpaceMaxCount : 200
EndPnts : 2.214
OctList : 3.33
SingleList : 1 : 4.128

Ver0.0.3
** start003 ** : 4100本
InSpaceMaxCount : 200
EndPnts : 1.794
OctList : 3.103
SingleList : 1 : 3.909

Ver0.0.2
** start ** : 8200本
OctLv : 7  InSpaceMaxCount : 200
EndPnts : 3.836
OctList : 7.519
SingleList : 2 : 10.692

Ver0.0.3
** start003 ** : 8200本
InSpaceMaxCount : 200
EndPnts : 4.359
OctList : 8.764
SingleList : 2 : 11.909

Ver0.0.2
** start ** : 12300本
OctLv : 7  InSpaceMaxCount : 200
EndPnts : 6.409
OctList : 14.139
SingleList : 3 : 21.447

Ver0.0.3
** start003 ** : 12300本
InSpaceMaxCount : 200
EndPnts : 6.472
OctList : 15.635
SingleList : 3 : 22.646

逆に1秒程遅くなっちゃいましたね･･･。

もう少し効率を上げられそうなｱｲﾃﾞｱは幾つかあるのですが、単独な3D曲線を
見つけ出すﾏｸﾛ自体が、恐らくご相談頂いた方以外には需要がほぼ無い気が
するので、この辺で一旦終了にします。

8分木空間分割ﾗｲﾌﾞﾗﾘについては結構悩みながら作ったので、他のﾏｸﾛでも
利用可能なように修正し、こちらのﾏｸﾛ達にも実装すれば処理時間の短縮が
出来そうな気がしています。こちら何かで。
ｻｰﾌｪｽの色をﾎﾞﾃﾞｨに反映する - C#ATIA
二つのﾎﾞﾃﾞｨ/形状ｾｯﾄを比較して、差分を抽出する - C#ATIA
そのうちに…。

ｺﾚｸｼｮﾝのﾊﾞﾌﾞﾙｿｰﾄをお借りする

前回が今年の最後にするつもりだったのですが、気になる事があったので
これを最後にします。

VBAでｺﾚｸｼｮﾝや配列をｿｰﾄしたい場合基本的には無く、単調なものであれば
DotNetのArrayListを利用する事も考えるのですが、ちょっと単調ではない
条件(ﾌﾟﾛﾊﾟﾃｨで　とか)でｿｰﾄする場合、ｺﾞﾘｺﾞﾘ書くしか方法がやはり無いん
ですね。

非常に参考になったのがこちらのt-hom'sさんのｻｲﾄ

http://thom.hateblo.jp/entry/2015/11/29/212934

これを流用させてもらおうと思ったのですが、案の定CATIAのApplicationｸﾗｽ
にはRunﾒｿｯﾄﾞ無いんですよ･･･。

好み的な部分も有り、こんな感じにしてみました。
(ﾃﾞｰﾀが整数って言うのが、無意味に近いのですが)

'vba　Collection_BubbleSort_Test

Sub Collection_BubbleSort_Test()
    Const Min = 1& '最小値
    Const Max = 100& '最大値
    Const Count = 10& '数
    
    Dim lst As Collection: Set lst = New Collection
    Dim i&
    For i = 0 To Count
        lst.Add Fix((Max - Min + 1) * Rnd + 1)
    Next
    
    Debug.Print "-- Before --"
    Call Dump(lst)
    Debug.Print "-- Aafter --"
    Call B_Sort_List(lst)
    Call Dump(lst)
    Debug.Print "-- End --"
End Sub

'ｿｰﾄ条件 True-降順 False-昇順
Private Function Comparison(ByVal A As Variant, ByVal B As Variant) As Boolean
    Comparison = A > B
End Function

'BubbleSort_Collection
Private Sub B_Sort_List(ByRef List As Collection)
    Dim i&, j&
    For i = 1 To List.Count
        For j = List.Count To i Step -1
            If Comparison(List(i), List(j)) Then
                    CollectionSwap List, i, j
            End If
        Next j
    Next i
End Sub

'Collection用ｽﾜｯﾌﾟ
Private Sub CollectionSwap(ByRef List As Collection, ByVal Idx1&, ByVal Idx2&)
    Dim Item1, Item2
    With List
        If IsObject(.Item(Idx1)) Then
            Set Item1 = .Item(Idx1)
            Set Item2 = .Item(Idx2)
        Else
            Let Item1 = .Item(Idx1)
            Let Item2 = .Item(Idx2)
        End If
        .Add Item1, After:=Idx2
        .Remove Idx2
        .Add Item2, After:=Idx1
        .Remove Idx1
    End With
End Sub

'Debug
Private Sub Dump(ByVal lst)
    Dim l
    For Each l In lst
        Debug.Print l
    Next
End Sub

キモは、Comparison関数です。ﾊﾞﾌﾞﾙｿｰﾄ内の IF文 で、True か False を返す
関数さえ作ってしまえば良いのかな?　と思いました。
上記では昇順ですが、降順にしたい場合はComparison関数内を

    Comparison = A < B

に、するだけです。

但し、t-hom'sさんのｻﾝﾌﾟﾙほど汎用性が無いんです。複数の異なる条件でｿｰﾄ
したい場合は、NGなんです･･･。

それを考えると Ariawase の Funcｸﾗｽを利用し

Private Sub B_Sort_List(ByRef List As Collection, ByVal fun As Func)
　・・・

の様にし、ｿｰﾄ条件用の関数さえ用意しておけば、異なる条件化で
ｿｰﾄ出来る様になるんじゃないのかな?　と思ったりしてます。(未確認です)
悩む･･･。　良いお年を。

単独な3D曲線の取得する2

こちらの続きです。
単独な3D曲線の取得する1 - C#ATIA

前回は総当りの非常に効率の悪い物でしたが、線形8分木による
ｽｸﾘｰﾆﾝｸﾞ処理を利用し単独な3D曲線を検索するように修正しました。

まずはCATMainを持つ標準ﾓｼﾞｭｰﾙ "FindSingle3DCurve.bas" です。

'vba FindSingle3DCurve Ver0.0.2  using-'KCL0.09' 'FindSingle3DCurve_Octree0.0.1'
'単独な3D曲線の取得

Option Explicit

'*** 設定 ***
Private Const Tolerance = 0.001                         '端点一致ﾄﾚﾗﾝｽ
Private Const OctLv = 2&                                '8分木分割ﾚﾍﾞﾙ(最大7)
Private Const InSpaceMaxCount = 200&                    '同一空間内最大数(目安)
Private Const SingleCurveFooter = " - SingleCurve!!"    '単独曲線名に追記する文字
'************

Sub CATMain()
    'ﾄﾞｷｭﾒﾝﾄのﾁｪｯｸ
    If Not CanExecute(Array("PartDocument", "ProductDocument")) Then Exit Sub
    Dim Msg$
    
    '形状ｾｯﾄの選択
    Msg = "3D曲線の入っている形状ｾｯﾄを指定してください : ESCｷｰ 終了"
    Dim Hb As HybridBody
    Set Hb = KCL.SelectItem(Msg, "HybridBody")
    If KCL.IsNothing(Hb) Then Exit Sub
    
    '線の取得
    Dim Lines As Collection: Set Lines = GetLines(Hb)
    If KCL.IsNothing(Lines) Then
        Msg = "指定した' " + Hb.Name + " 'には、3D曲線がありませんでした!"
        MsgBox Msg: Exit Sub
    End If
    
    '警告
    Msg = CStr(Lines.Count) + "本あります。実行しますか?"
    
    If MsgBox(Msg, vbYesNo + vbInformation) = vbNo Then Exit Sub
    Debug.Print "** start ** : " + CStr(Lines.Count) + "本" 'debug
    Debug.Print "OctLv : " + CStr(OctLv) + "  InSpaceMaxCount : " + CStr(InSpaceMaxCount) 'debug
    KCL.SW_Start 'debug
    
    '端点座標取得
    Dim EndPnts As Collection: Set EndPnts = GetEndPntCoordinates(Lines)
    Debug.Print "EndPnts : " + CStr(KCL.SW_GetTime) 'debug
    
    '組み合わせﾘｽﾄ取得
    Dim OctList As Collection
    If FindSingle3DCurve_Octree.Init(OctLv, Tolerance, InSpaceMaxCount) Then
        Set OctList = FindSingle3DCurve_Octree.GetLinerOctreeList(EndPnts)
    End If
    
    If KCL.IsNothing(OctList) Then
        Msg = "ｴﾗｰです。ﾏｸﾛ内の設定値を修正して下さい(特に 'OctLv' )"
        MsgBox Msg: Exit Sub
    End If
    Debug.Print "OctList : " + CStr(KCL.SW_GetTime) 'debug
    
    '単独曲線ｲﾝﾃﾞｯｸｽ取得
    Dim SingleList As Collection: Set SingleList = GetSingleIdx(EndPnts, OctList)
    Debug.Print "SingleList : " + CStr(SingleList.Count) + " : " + CStr(KCL.SW_GetTime)
    If SingleList.Count < 1 Then
        Msg = "指定した' " + Hb.Name + " 'には、単独な3D曲線が見つかりませんでした!" + vbNewLine + _
              "(ﾄﾚﾗﾝｽ " & Tolerance & " mm)"
        MsgBox Msg: Exit Sub
    End If
    
    '確認
    Msg = SingleList.Count & "本の単独な3D曲線が見つかりました!!" + vbNewLine + _
          "名前の最後に ' " + SingleCurveFooter + " ' を追記しますか?"
    If MsgBox(Msg, vbYesNo + vbInformation) = vbNo Then GoTo MacroEnd
    
    'ﾘﾈｰﾑ
    Call AppendName(Lines, SingleList)
    
MacroEnd:
    'Call Dump(ByVal SingleList, ByVal Lines, ByVal EndPnts) 'debug
    MsgBox "終了!"
End Sub

'単独曲線のﾘﾈｰﾑ
Private Sub AppendName(ByVal Lines As Collection, ByVal IdxList As Collection)
    Dim i
    For Each i In IdxList
        With Lines.Item(i)
            .Name = .Name & SingleCurveFooter
        End With
    Next
End Sub

'初期化済み配列生成 - ｵﾌﾞｼﾞｪｸﾄNG
Private Function InitRangeAry(ByVal Count&, ByVal Value As Variant)
    Dim ary() As Variant: ReDim ary(Count)
    Dim i&
    For i = 0 To Count
        ary(i) = Value
    Next
    InitRangeAry = ary
End Function

'単独曲線ｲﾝﾃﾞｯｸｽ取得
Private Function GetSingleIdx(ByVal EndPnts As Collection, ByVal OctList As Collection) As Collection
    Dim i&, j&, i_crv_idx, j_crv_idx
    Dim CombAry As Variant: CombAry = InitRangeAry(EndPnts.Count * 0.5, False)
    Dim Space As Collection
    For Each Space In OctList
        For i = 1 To Space.Count
            i_crv_idx = EndPnts.Item(Space.Item(i))(3)
            For j = i + 1 To Space.Count
                j_crv_idx = EndPnts.Item(Space.Item(j))(3)
                If i_crv_idx = j_crv_idx Then GoTo Continue
                If CombAry(i_crv_idx) And CombAry(j_crv_idx) Then GoTo Continue
                If IsPosEqual(EndPnts.Item(Space.Item(i)), EndPnts.Item(Space.Item(j))) Then
                    CombAry(i_crv_idx) = True
                    CombAry(j_crv_idx) = True
                End If
Continue:
            Next
        Next
    Next
    
    Dim SingleIdxList As Collection: Set SingleIdxList = New Collection
    For i = 1 To UBound(CombAry)
        If CombAry(i) = False Then SingleIdxList.add i
    Next
    Set GetSingleIdx = SingleIdxList
End Function

'座標値一致
Private Function IsPosEqual(ByVal p1 As Variant, ByVal p2 As Variant) As Boolean
    IsPosEqual = False
    Dim i&
    For i = 0 To 2
        If Abs(p2(i) - p1(i)) > Tolerance Then Exit Function
    Next
    If Sqr(Abs((p2(0) - p1(0)) * (p2(0) - p1(0)) + _
               (p2(1) - p1(1)) * (p2(1) - p1(1)) + _
               (p2(2) - p1(2)) * (p2(2) - p1(2)))) < Tolerance Then
        IsPosEqual = True
    End If
End Function

'端点取得
Private Function GetEndPntCoordinates(ByRef Lines As Collection) As Collection
    Dim Doc As PartDocument: Set Doc = KCL.GetParent_Of_T(Lines.Item(1), "PartDocument")
    Dim SPAWB As Workbench: Set SPAWB = Doc.GetWorkbench("SPAWorkbench")
    Dim EndPnt As Collection: Set EndPnt = New Collection
    Dim Pos(8) As Variant
    Dim i&
    For i = 1 To Lines.Count
        Call SPAWB.GetMeasurable(Lines.Item(i)).GetPointsOnCurve(Pos)
        EndPnt.add KCL.JoinAry(KCL.GetRangeAry(Pos, 0, 2), Array(i))
        EndPnt.add KCL.JoinAry(KCL.GetRangeAry(Pos, 6, 8), Array(i))
    Next
    Set GetEndPntCoordinates = EndPnt
End Function

'線取得
Private Function GetLines(ByVal Hb As HybridBody) As Collection
    Dim Hss As HybridShapes: Set Hss = Hb.HybridShapes
    If Hss.Count < 1 Then Exit Function
    
    Dim Pt As Part: Set Pt = KCL.GetParent_Of_T(Hb, "PartDocument").Part
    Dim Fact As HybridShapeFactory: Set Fact = Pt.HybridShapeFactory
    Dim Lines As Collection: Set Lines = New Collection
    Dim Hs As HybridShape
    For Each Hs In Hss
        Select Case Fact.GetGeometricalFeatureType(Hs)
            Case 2, 3, 4
                Lines.add Hs
        End Select
    Next
    Set GetLines = Lines
End Function

'debug用
Private Sub Dump(ByVal SingleList, ByVal Lines, ByVal EndPnts)
    Dim ss
    For Each ss In SingleList
        Debug.Print ss & ":" & Lines.Item(EndPnts.Item(ss)(3)).Name
    Next
End Sub

'ver0.0.1 - 総当りのﾀﾞﾒなやつ
'ver0.0.2 - 線形8分木を実装

設定を行う為の定数部分の説明です。
・"Tolerance"
　端点一致を判断するﾄﾚﾗﾝｽです。一般的に考え0.001以下は無意味な気がしています。
　又、ﾃﾞｰﾀの状態次第ですが、空間分割する際にも影響がある為あまり
　大きな値(100.0とか)では、処理時間が長くなる可能性があります。

・"OctLv"
　空間分割する際の分割ﾚﾍﾞﾙです。ﾃｽﾄを行った感じでは大きくした方が処理時間が短く
　なる傾向にあります。　最大7です。

・"InSpaceMaxCount"
　空間分割した際、同一空間内の要素数の最大数ですが、あくまで目安です。
　こちらは値は、大小どちらにした方が効率良くなるのか正直把握しておりません。
　
・"SingleCurveFooter"
　前回同様、Hitした単独曲線名に追記する文字です。

続いて、組み合わせ効率を改善するための標準ﾓｼﾞｭｰﾙ "FindSingle3DCurve_Octree.bas" です。
あまりに複雑な上、他にも流用できそうな気がしたので、単独なﾓｼﾞｭｰﾙとしました。

'vba FindSingle3DCurve_Octree Ver0.0.1  using-'KCL0.08'
'ﾓｰﾄﾝ順序を利用した8分木空間分割ﾗｲﾌﾞﾗﾘ
'FindSingle3DCurveﾏｸﾛ専用です

Option Explicit

Private Const CLINER8TREEMANAGER_MAXLEVEL = 7   '有効空間分割最大ﾚﾍﾞﾙ

Private m_Level&                                '分割ﾚﾍﾞﾙ
Private m_Tolerance#                            '一致ﾄﾚﾗﾝｽ
Private m_MaxCount&                             '同一空間内最大数(目安)
Private m_AxisCount&                            '空間分割時の各軸の最大数
Private m_MinPos                                '空間最小座標
Private m_Unit                                  '空間単位ｻｲｽﾞ
Private m_ToleranceRatio                        '空間単位ｻｲｽﾞに対してのﾄﾚﾗﾝｽ比率
Private m_CellCount&()                          'ﾚﾍﾞﾙ毎の空間数

'*** Octree ***
'http://marupeke296.com/COL_3D_No15_Octree.html

'線形8分木準備
''' @param :Level-long-分割ﾚﾍﾞﾙ
''' @param :Tolerance-Double-一致ﾄﾚﾗﾝｽ
''' @param :MaxCount-long-同一空間内最大数(目安)
''' @return:Boolean
Function Init(ByVal Level&, ByVal Tolerance#, ByVal MaxCount&) As Boolean
    Init = False
    If Level > CLINER8TREEMANAGER_MAXLEVEL Then Exit Function
    If Tolerance <= 0 Then Exit Function
    
    m_Level = Level
    m_Tolerance = Tolerance
    m_MaxCount = MaxCount
    m_ToleranceRatio = InitRangeAry(2, 0)
    m_AxisCount = sl(1, m_Level)
    
    ReDim m_CellCount(CLINER8TREEMANAGER_MAXLEVEL + 1)
    m_CellCount(0) = 1
    Dim i&
    For i = 1 To UBound(m_CellCount)
        m_CellCount(i) = m_CellCount(i - 1) * 8
    Next
    Init = True
End Function

'線形8分木ﾘｽﾄ取得
''' @param :Pnts-Collection(array(Double))-座標値郡
''' @return:Collection(Collection(long))-空間別点Idx郡
Function GetLinerOctreeList(ByVal Pnts As Collection)
    '座標値郡のIdxList作成
    Dim PntIdxList As Collection: Set PntIdxList = InitRangeList(Pnts.Count)
    
    '空間内最大数以下の場合、そのまま返す
    Dim DecidedList As Collection: Set DecidedList = New Collection
    If Pnts.Count < m_MaxCount Then
        Call DecidedList.add(PntIdxList)
        GoTo FuncEnd
    End If
    
    'IdxListをﾙｰﾄ空間として登録
    Dim CheckList As Collection: Set CheckList = New Collection
    Call CheckList.add(PntIdxList)
    
    'ﾚﾍﾞﾙの応じたIdx用配列
    Dim SpeceEnum, TempList As Collection, ReCheckList As Collection
    Dim Space, TempSpace, SpaAry, i&, Idx
    Do
        Set ReCheckList = New Collection
        '配置
        For Each Space In CheckList
            SpeceEnum = InitRangeAry(m_CellCount(m_Level), -1)
            Set TempList = New Collection
            
            If Not SetSpaceInfo(Pnts, Space) Then '空間が小さすぎる
                Call DecidedList.add(Space)
                GoTo Continue
            End If
            
            For Each Idx In Space
                SpaAry = GetMortonNum(Pnts(Idx))
                For i = 0 To UBound(SpaAry)
                    If SpeceEnum(SpaAry(i)) < 0 Then
                        Call TempList.add(InitSpace())
                        SpeceEnum(SpaAry(i)) = TempList.Count
                    End If
                    Call TempList.Item(SpeceEnum(SpaAry(i))).add(Idx)
                Next
            Next
            
            '空間毎の数が多いものを再配置
            For Each TempSpace In TempList
                Select Case True
                    Case TempSpace.Count < 2
                        '空間に1個しかないものは検証しない
                        'Call DecidedList.add(TempSpace)
                    Case TempSpace.Count > m_MaxCount
                        Call ReCheckList.add(TempSpace)
                    Case Else
                        Call DecidedList.add(TempSpace)
                End Select
            Next
Continue:
        Next
        If ReCheckList.Count < 1 Then Exit Do
        Set CheckList = ReCheckList
    Loop
    
FuncEnd:
    Set GetLinerOctreeList = DecidedList
    Set PntIdxList = Nothing
    Set DecidedList = Nothing
    Set CheckList = Nothing
    Set ReCheckList = Nothing
    Set TempList = Nothing
    ReDim SpeceEnum(0)
End Function

'空のｺﾚｸｼｮﾝ
Private Function InitSpace() As Collection
    Set InitSpace = New Collection
End Function

'空間情報設定
''' @param :Pnts-Collection(array(Double))-座標値郡
''' @param :Idxs-Collection(long)-座標値郡Idx
''' @return:Boolean
Private Function SetSpaceInfo(ByVal Pnts, ByVal Idxs As Collection) As Boolean
    SetSpaceInfo = False
    Dim SpSize: SpSize = GetSpaceSize_Idx(Pnts, Idxs)
    m_MinPos = AryAdd(SpSize(0), m_Tolerance * -1)
    Dim w: w = ArySub(AryAdd(SpSize(1), m_Tolerance), m_MinPos)
    m_Unit = AryDiv(w, m_AxisCount)
    
    Dim MinUnit#: MinUnit = m_Tolerance * 3
    Dim i&
    For i = 0 To 2
         If m_Unit(i) < MinUnit Then Exit Function
    Next
    
    For i = 0 To 2
        m_ToleranceRatio(i) = m_ToleranceRatio(i) / m_Unit(i)
    Next
    SetSpaceInfo = True
End Function

'座標→最小ﾚﾍﾞﾙ空間番号郡取得
''' @param :Pos-array(Double)-座標値
''' @return:array(Long)
Private Function GetMortonNum(ByVal Pos As Variant) As Variant
    Dim Ratio#(2), Inte&(2), Dec#(2)
    Dim i&
    For i = 0 To 2
        Ratio(i) = (Pos(i) - m_MinPos(i)) / m_Unit(i)
        Inte(i) = Fix(Ratio(i))
        Dec(i) = Ratio(i) - Inte(i)
    Next
    
    Dim Axis(2) As Variant
    Dim AxisNums As Collection
    For i = 0 To 2
        Set AxisNums = New Collection
        AxisNums.add Inte(i)
        If Dec(i) <= m_ToleranceRatio(i) And Inte(i) > 0 Then AxisNums.add Inte(i) - 1
        If 1# - Dec(i) <= m_ToleranceRatio(i) And Inte(i) < m_AxisCount - 1 Then AxisNums.add Inte(i) + 1
        Set Axis(i) = AxisNums
    Next
    
    Dim x, y, z
    Dim SpaNo(): ReDim SpaNo(Axis(0).Count * Axis(1).Count * Axis(2).Count)
    Dim cnt&: cnt = -1
    For Each x In Axis(0)
        For Each y In Axis(1)
            For Each z In Axis(2)
                cnt = cnt + 1
                SpaNo(cnt) = Get3DMortonNumber(x, y, z)
                If SpaNo(cnt) < 0 And SpaNo(cnt) >= m_CellCount(m_Level) Then cnt = cnt - 1
    Next: Next: Next
    ReDim Preserve SpaNo(cnt)
    GetMortonNum = SpaNo
End Function

'ビット分割関数
''' @param :n-long
''' @return:long
Private Function BitSeparateFor3D(ByVal n&) As Long
    Dim s As Long: s = n
    s = (s Or sl(s, 8)) And &HF00F
    s = (s Or sl(s, 4)) And &HC30C3
    s = (s Or sl(s, 2)) And &H249249
    BitSeparateFor3D = s
End Function

'8分木モートン順序算出関数
''' @param :x-long
''' @param :y-long
''' @param :z-long
''' @return:long
Private Function Get3DMortonNumber(ByVal x&, ByVal y&, ByVal z&) As Long
   Get3DMortonNumber = BitSeparateFor3D(x) Or _
                       sl(BitSeparateFor3D(y), 1) Or _
                       sl(BitSeparateFor3D(z), 2)
End Function

'座標値郡から空間ｻｲｽﾞ取得
''' @param :EndPnts-Collection(array(Double))-座標値郡
''' @return:array(array(Double))-0:最小値 1:最大値
Private Function GetSpaceSize_Idx(ByVal EndPnts As Collection, ByVal IdxList As Collection) As Variant
    Dim Min: Min = InitRangeAry(2, 1.79769313486231E+308)
    Dim Max: Max = InitRangeAry(2, -1.79769313486231E+308)
    
    Dim Idx, i&
    For Each Idx In IdxList
        For i = 0 To 2
            If Min(i) > EndPnts.Item(Idx)(i) Then Min(i) = EndPnts.Item(Idx)(i)
            If Max(i) < EndPnts.Item(Idx)(i) Then Max(i) = EndPnts.Item(Idx)(i)
        Next
    Next
    GetSpaceSize_Idx = Array(Min, Max)
End Function

'*** BitShift ***
'http://www.geocities.co.jp/SiliconValley/4334/unibon/asp/bitshift2.html
' 左シフト
Private Function sl(ByVal x&, ByVal n&) As Long
    If n = 0 Then
        sl = x
    Else
        Dim k: k = CLng(2 ^ (32 - n - 1))
        Dim d: d = x And (k - 1)
        Dim c: c = d * CLng(2 ^ n)
        If x And k Then c = c Or &H80000000
        sl = c
    End If
End Function

'*** VBA不足関数 ***
'配列同士の引き算-細かいチェック無し
Private Function ArySub(ByVal a, ByVal b) As Variant
    ArySub = Array(a(0) - b(0), a(1) - b(1), a(2) - b(2))
End Function

'配列と実数の足し算-細かいチェック無し
Private Function AryAdd(ByVal a, ByVal b#) As Variant
    AryAdd = Array(a(0) + b, a(1) + b, a(2) + b)
End Function

'配列と実数の割り算-細かいチェック無し
Private Function AryDiv(ByVal a, ByVal b#) As Variant
    AryDiv = Array(a(0) / b, a(1) / b, a(2) / b)
End Function

'初期化済み配列生成 - ｵﾌﾞｼﾞｪｸﾄNG
Private Function InitRangeAry(ByVal Count&, ByVal Value As Variant)
    Dim ary() As Variant: ReDim ary(Count)
    Dim i&
    For i = 0 To Count
        ary(i) = Value
    Next
    InitRangeAry = ary
End Function

'初期化済みｺﾚｸｼｮﾝ生成
Private Function InitRangeList(ByVal Count&) As Collection
    Dim List As Collection: Set List = New Collection
    Dim i&
    For i = 1 To Count
        List.add i
    Next
    Set InitRangeList = List
End Function

これは単独では利用しません。かなり解りにくいのですが、一応説明です。
Init関数に、分割ﾚﾍﾞﾙ,・ﾄﾚﾗﾝｽ,・空間内最大数 を渡して初期化し
GetLinerOctreeList関数に座標値の入ったｺﾚｸｼｮﾝを渡すことで、
組み合わせを行うための、座標値のｲﾝﾃﾞｯｸｽが空間毎に入ったｺﾚｸｼｮﾝを返します。
引数-Collection(Array(x,y,z-Double)) の形で配列は3次元以上でもOKです。
戻り値-Collection(Collection(Long))です。返すのはあくまでも列挙する為の
ｲﾝﾃﾞｯｸｽで、座標値ではありません。

後は、標準ﾓｼﾞｭｰﾙ KCL が必要になります。

こんな感じで3つが同一ﾌﾟﾛｼﾞｪｸﾄに入っている必要があります。

これを、前回(Ver0.0.1)の物と比較しながら実際に処理した結果です。

Ver0.0.1
** start ** : 4100本
EndPnts : 1.73
SingleList : 1 : 2771.944

Ver0.0.2
** start ** : 4100本
OctLv : 2  InSpaceMaxCount : 200
EndPnts : 1.771
OctList : 3.607
SingleList : 1 : 9.34

Ver0.0.2
** start ** : 4100本
OctLv : 7  InSpaceMaxCount : 200
EndPnts : 2.214
OctList : 3.33
SingleList : 1 : 4.128

SingleList行の最後の値が、ﾄｰﾀﾙの処理時間です。
4099本→4100本に1本増やしたのは、使用していたｻﾝﾌﾟﾙﾃﾞｰﾀに単独な
3D曲線が無かった為です。

4100本を処理した際、Ver0.0.1では45分以上の処理時間が必要だったのですが
Ver0.0.2では、分割ﾚﾍﾞﾙ2で9秒　分割ﾚﾍﾞﾙ7では4秒程まで短縮しました。
線形8分木の効果が絶大です。　これであれば十分耐えられる処理時間かと
思います。
但し、"FindSingle3DCurve.bas" 側ではもう少し効率化出来そうなｱｲﾃﾞｱが
幾つかあるので、今後の課題とします。

正直に書くと、線形8分木に取り組んだのは初めてでしたが、線形4分木に
ついては過去にDrawの重複線削除ﾏｸﾛで取り組んだことがあります。
ですが、線形4分木の衝突ﾘｽﾄを作成する為の空間の巡回方法がどうしても
理解できず、情けない話断念した経緯があり、今回も実装までたどり着けるか
どうか当初は解りませんでした。
(このﾌﾞﾛｸﾞを始めた理由の一つは、重複線削除ﾏｸﾛをC#化し線形4分木を
実装することでした。 …全くやっていませんが)

線形8分木の実装は年内では間に合わないと思っていたのですが、こちらの
思い付きの効果が大きく
端点一致を探る、組み合わせテスト10 - C#ATIA
ｷﾞﾘｷﾞﾘですが、年内には出来上がったつもりです。(ﾊﾞｸﾞがあるかも…)

今年はこれでおしまいです。昨年も書きましたが、こんなにﾌﾞﾛｸﾞを続けると
思っていませんでしたよ。

不覚にも、KCLのご利用をお考えの方へ

非常に個人的に作成している CATVBA用ﾗｲﾌﾞﾗﾘ なのですが、ｻﾝﾌﾟﾙとして
掲載しているﾏｸﾛを利用する際には必要となる為、今更ながらご説明です。

例えばこちらの場合ですが、
単独な3D曲線の取得する1 - C#ATIA

ｻﾝﾌﾟﾙｺｰﾄﾞのﾍｯﾀﾞ部分を確認してください。
(忘れていなければ) ﾀｸﾞが有りますので、ｸﾘｯｸして下さい。

過去ﾛｸﾞが表示されます。一番古いﾛｸﾞにｺｰﾄﾞが掲載されています。
基本的にKCLに関しては、ここを書き換えることにしております。

非常に個人的なCATVBA用ﾗｲﾌﾞﾗﾘ - C#ATIA
こちらが最新のものとなります。
GitHub - kantoku-code/KCL: Library for personal CATVBA (CATIA macro)

大きな機能追加についてはｱﾅｳﾝｽしようかとは考えておりますが、
利用に関して影響の無い修正ﾚﾍﾞﾙのﾊﾞｰｼﾞｮﾝUpについては
特にｱﾅｳﾝｽしません。(大半が当方のﾎﾟｶﾐｽの為です･･･)
その為、ｻﾝﾌﾟﾙｺｰﾄﾞに記載されているﾊﾞｰｼﾞｮﾝにも注意しておいて下さい。

先程のKCLのｺｰﾄﾞをｺﾋﾟｰしてもらい、新規の標準ﾓｼﾞｭｰﾙ作成しﾍﾞﾀｯと
ﾍﾟｰｽﾄして頂き、ﾓｼﾞｭｰﾙ名を "KCL" としておいて下さい。
ｻﾝﾌﾟﾙｺｰﾄﾞについても同様に、新規の標準ﾓｼﾞｭｰﾙ作成しｺﾋﾟﾍﾟして下さい。
(ﾓｼﾞｭｰﾙ名はご自由に。あくまでも同一ﾌﾟﾛｼﾞｪｸﾄで行ってください。)

･･･最近ずっと悩んでいるのですが、' Ariawase ' には、非常に魅力的な
配列操作の関数郡が有り、KCLのクズ配列操作関数を全て削除して
Ariawaseを併用したいような気持ちが･･･。
過去にFuncｸﾗｽがCATIAでは上手く動作しなかったのですが、複数の戻り値が
返せないVBAとしては、TupleとかArrZipが欲しい。かなり迷う。

端点一致を探る、組み合わせテスト11

こちらの続きです。
端点一致を探る、組み合わせテスト10 - C#ATIA

前回の物を実装する為に書き換えました。

'vba
'ﾓｰﾄﾝ順序を利用した8分木空間分割の為のｸｽﾞｺｰﾄﾞ4
Option Explicit

Private Const CLINER8TREEMANAGER_MAXLEVEL = 7   '有効空間分割最大ﾚﾍﾞﾙ

Private m_Level&                                '分割ﾚﾍﾞﾙ
Private m_Tolerance#                            '一致ﾄﾚﾗﾝｽ
Private m_MaxCount&                             '同一空間内最大数(目安)
Private m_AxisCount&                            '空間分割時の各軸の最大数
Private m_MinPos                                '空間最小座標
Private m_Unit                                  '空間単位ｻｲｽﾞ
Private m_ToleranceRatio                        '空間単位ｻｲｽﾞに対してのﾄﾚﾗﾝｽ比率
Private m_CellCount&()                          'ﾚﾍﾞﾙ毎の空間数

Sub CATMain()
    'ﾃｽﾄ座標
    Dim Poss As Collection: Set Poss = New Collection
    With Poss
        .add Array(-40#, -40#, -40#)
        .add Array(40#, 40#, 40#)
        .add Array(0, 0, 0)
        .add Array(0.1, 0, 0)
        .add Array(0.1, 0.001, 0)
        .add Array(-39, -39, -39)
        .add Array(-39, -39.1, -39)
        .add Array(-39, -39, -39.1)
    End With
    
    Dim lv&: lv = 2                  '分割ﾚﾍﾞﾙ
    Dim Tol#: Tol = 0.001            '端点一致ﾄﾚﾗﾝｽ
    Dim MaxCount&: MaxCount = 3      '空間内最大数
    
    If Not Init(lv, Tol, MaxCount) Then Exit Sub
    
    '分割した空間毎に組み合わせる端点のｲﾝﾃﾞｯｸｽを取得
    Dim List As Collection: Set List = GetLinerOctreeList(Poss)
    
    Dim msg$
    msg = "分割ﾚﾍﾞﾙ : " & lv & vbNewLine + _
          "端点一致ﾄﾚﾗﾝｽ : " & Tol & vbNewLine + _
          "空間内最大数 : " & MaxCount & vbNewLine + _
          "空間数 : " & List.Count & vbNewLine + _
          "以下、組み合わせるﾘｽﾄです" & vbNewLine
    Dim msg2$, tmp: msg2 = ""
    For Each tmp In List
        msg2 = msg2 & List2Str(tmp) & vbNewLine
    Next
    MsgBox msg & msg2
End Sub

'ﾃｽﾄ用
Private Function List2Str(ByVal List) As String
    Dim tmp, msg$: msg = ""
    For Each tmp In List
        msg = msg & tmp & ","
    Next
    List2Str = Mid(msg, 1, Len(msg) - 1)
End Function


'*** Octree ***
'http://marupeke296.com/COL_3D_No15_Octree.html

'線形8分木準備
''' @param :Level-long-分割ﾚﾍﾞﾙ
''' @param :Tolerance-Double-一致ﾄﾚﾗﾝｽ
''' @param :MaxCount-long-同一空間内最大数(目安)
''' @return:Boolean
Function Init(ByVal Level&, ByVal Tolerance#, ByVal MaxCount&) As Boolean
    Init = False
    If Level >= CLINER8TREEMANAGER_MAXLEVEL Then Exit Function
    If Tolerance <= 0 Then Exit Function
    
    m_Level = Level
    m_Tolerance = Tolerance
    m_MaxCount = MaxCount
    m_ToleranceRatio = InitRangeAry(2, 0)
    m_AxisCount = sl(1, m_Level)
    
    ReDim m_CellCount(CLINER8TREEMANAGER_MAXLEVEL + 1)
    m_CellCount(0) = 1
    Dim i&
    For i = 1 To UBound(m_CellCount)
        m_CellCount(i) = m_CellCount(i - 1) * 8
    Next
    Init = True
End Function

'線形8分木ﾘｽﾄ取得
''' @param :Pnts-Collection(array(Double))-座標値郡
''' @return:Collection(Collection(long))-空間別点Idx郡
Function GetLinerOctreeList(ByVal Pnts As Collection)
    '空間内最大数以下の場合、そのまま返す
    Dim DecidedList As Collection: Set DecidedList = New Collection
    If Pnts.Count < m_MaxCount Then
        Call DecidedList.add(Pnts)
        GoTo FuncEnd
    End If
    
    '座標値郡のIdxList作成
    Dim PntIdxList As Collection: Set PntIdxList = InitRangeList(Pnts.Count)
    
    'IdxListをﾙｰﾄ空間として登録
    Dim CheckList As Collection: Set CheckList = New Collection
    Call CheckList.add(PntIdxList)
    
    'ﾚﾍﾞﾙの応じたIdx用配列
    Dim SpeceEnum, TempList As Collection, ReCheckList As Collection
    Dim Space, TempSpace, SpaAry, i&, Idx
    Do
        Set ReCheckList = New Collection
        '配置
        For Each Space In CheckList
            SpeceEnum = InitRangeAry(m_CellCount(m_Level), -1)
            Set TempList = New Collection
            
            If Not SetSpaceInfo(Pnts, Space) Then '空間が小さすぎる
                Call DecidedList.add(Space)
                GoTo Continue
            End If
            
            For Each Idx In Space
                SpaAry = GetMortonNum(Pnts(Idx))
                For i = 0 To UBound(SpaAry)
                    If SpeceEnum(SpaAry(i)) < 0 Then
                        Call TempList.add(InitSpace())
                        SpeceEnum(SpaAry(i)) = TempList.Count
                    End If
                    Call TempList.Item(SpeceEnum(SpaAry(i))).add(Idx)
                Next
            Next
            
            '空間毎の数が多いものを再配置
            For Each TempSpace In TempList
                Select Case True
                    Case TempSpace.Count < 2
                        '空間に1個しかないものは検証しない
                    Case TempSpace.Count > m_MaxCount
                        Call ReCheckList.add(TempSpace)
                    Case Else
                        Call DecidedList.add(TempSpace)
                End Select
            Next
Continue:
        Next
        If ReCheckList.Count < 1 Then Exit Do
        Set CheckList = ReCheckList
    Loop
    
FuncEnd:
    Set GetLinerOctreeList = DecidedList
End Function

'空のｺﾚｸｼｮﾝ
Private Function InitSpace() As Collection
    Set InitSpace = New Collection
End Function

'空間情報設定
''' @param :Pnts-Collection(array(Double))-座標値郡
''' @param :Idxs-Collection(long)-座標値郡Idx
''' @return:Boolean
Private Function SetSpaceInfo(ByVal Pnts, ByVal Idxs As Collection) As Boolean
    SetSpaceInfo = False
    Dim SpSize: SpSize = GetSpaceSize_Idx(Pnts, Idxs)
    m_MinPos = AryAdd(SpSize(0), m_Tolerance * -1)
    Dim w: w = ArySub(AryAdd(SpSize(1), m_Tolerance), m_MinPos)
    m_Unit = AryDiv(w, m_AxisCount)
    
    Dim MinUnit#: MinUnit = m_Tolerance * 3
    Dim i&
    For i = 0 To 2
         If m_Unit(i) < MinUnit Then Exit Function
    Next
    
    For i = 0 To 2
        m_ToleranceRatio(i) = m_ToleranceRatio(i) / m_Unit(i)
    Next
    SetSpaceInfo = True
End Function

'座標→最小ﾚﾍﾞﾙ空間番号郡取得
''' @param :Pos-array(Double)-座標値
''' @return:array(Long)
Private Function GetMortonNum(ByVal Pos As Variant) As Variant
    Dim Ratio#(2), Inte&(2), Dec#(2)
    Dim i&
    For i = 0 To 2
        Ratio(i) = (Pos(i) - m_MinPos(i)) / m_Unit(i)
        Inte(i) = Fix(Ratio(i))
        Dec(i) = Ratio(i) - Inte(i)
    Next
    
    Dim Axis(2) As Variant
    Dim AxisNums As Collection
    For i = 0 To 2
        Set AxisNums = New Collection
        AxisNums.add Inte(i)
        If Dec(i) <= m_ToleranceRatio(i) And Inte(i) > 0 Then AxisNums.add Inte(i) - 1
        If 1# - Dec(i) <= m_ToleranceRatio(i) And Inte(i) < m_AxisCount - 1 Then AxisNums.add Inte(i) + 1
        Set Axis(i) = AxisNums
    Next
    
    Dim x, y, z
    Dim SpaNo(): ReDim SpaNo(Axis(0).Count * Axis(1).Count * Axis(2).Count)
    Dim cnt&: cnt = -1
    For Each x In Axis(0)
        For Each y In Axis(1)
            For Each z In Axis(2)
                cnt = cnt + 1
                SpaNo(cnt) = Get3DMortonNumber(x, y, z)
                If SpaNo(cnt) < 0 And SpaNo(cnt) >= m_CellCount(m_Level) Then cnt = cnt - 1
    Next: Next: Next
    ReDim Preserve SpaNo(cnt)
    GetMortonNum = SpaNo
End Function

'ビット分割関数
''' @param :n-long
''' @return:long
Private Function BitSeparateFor3D(ByVal n&) As Long
    Dim s As Long: s = n
    s = (s Or sl(s, 8)) And &HF00F
    s = (s Or sl(s, 4)) And &HC30C3
    s = (s Or sl(s, 2)) And &H249249
    BitSeparateFor3D = s
End Function

'8分木モートン順序算出関数
''' @param :x-long
''' @param :y-long
''' @param :z-long
''' @return:long
Private Function Get3DMortonNumber(ByVal x&, ByVal y&, ByVal z&) As Long
   Get3DMortonNumber = BitSeparateFor3D(x) Or _
                       sl(BitSeparateFor3D(y), 1) Or _
                       sl(BitSeparateFor3D(z), 2)
End Function

'座標値郡から空間ｻｲｽﾞ取得
''' @param :EndPnts-Collection(array(Double))-座標値郡
''' @return:array(array(Double))-0:最小値 1:最大値
Private Function GetSpaceSize_Idx(ByVal EndPnts As Collection, ByVal IdxList As Collection) As Variant
    Dim Min: Min = InitRangeAry(2, 1.79769313486231E+308)
    Dim Max: Max = InitRangeAry(2, -1.79769313486231E+308)
    
    Dim Idx, i&
    For Each Idx In IdxList
        For i = 0 To 2
            If Min(i) > EndPnts.Item(Idx)(i) Then Min(i) = EndPnts.Item(Idx)(i)
            If Max(i) < EndPnts.Item(Idx)(i) Then Max(i) = EndPnts.Item(Idx)(i)
        Next
    Next
    GetSpaceSize_Idx = Array(Min, Max)
End Function

'*** BitShift ***
'http://www.geocities.co.jp/SiliconValley/4334/unibon/asp/bitshift2.html
' 左シフト
Private Function sl(ByVal x&, ByVal n&) As Long
    If n = 0 Then
        sl = x
    Else
        Dim k: k = CLng(2 ^ (32 - n - 1))
        Dim d: d = x And (k - 1)
        Dim c: c = d * CLng(2 ^ n)
        If x And k Then c = c Or &H80000000
        sl = c
    End If
End Function

'*** VBA不足関数 ***
'配列同士の引き算-細かいチェック無し
Private Function ArySub(ByVal a, ByVal b) As Variant
    ArySub = Array(a(0) - b(0), a(1) - b(1), a(2) - b(2))
End Function

'配列と実数の足し算-細かいチェック無し
Private Function AryAdd(ByVal a, ByVal b#) As Variant
    AryAdd = Array(a(0) + b, a(1) + b, a(2) + b)
End Function

'配列と実数の割り算-細かいチェック無し
Private Function AryDiv(ByVal a, ByVal b#) As Variant
    AryDiv = Array(a(0) / b, a(1) / b, a(2) / b)
End Function

'初期化済み配列生成 - ｵﾌﾞｼﾞｪｸﾄNG
Private Function InitRangeAry(ByVal Count&, ByVal Value As Variant)
    Dim ary() As Variant: ReDim ary(Count)
    Dim i&
        For i = 0 To Count
            ary(i) = Value
        Next
    InitRangeAry = ary
End Function

'初期化済みｺﾚｸｼｮﾝ生成
Private Function InitRangeList(ByVal Count&) As Collection
    Dim List As Collection: Set List = New Collection
    Dim i&
    For i = 1 To Count
        List.add i
    Next
    Set InitRangeList = List
End Function

Init関数で線形8分木を作成する為の準備を行い、GetLinerOctreeList関数に
座標値のｺﾚｸｼｮﾝを渡す事で、空間分割された点群のｲﾝﾃﾞｯｸｽを返します。

もう、GetLinerOctreeList関数の酷さ(3重ﾙｰﾌﾟ)と言ったら、言葉も出ない程
なのですが、他に思いつきませんでした･･･。

実際に実行した結果はこんな感じです。

他の方には意味不明かと思いますが、実はｵﾘｼﾞﾅﾙより効率が悪い部分が有ります。
両方のｱﾙｺﾞﾘｽﾞﾑを上手く利用できると効率良く処理できるのですが･･･。

ﾃｽﾄ段階ですが、こちらのｺｰﾄﾞに今回の物を実装したのですが

単独な3D曲線の取得する1 - C#ATIA

劇的に速くなりました。　線形8分木の素晴らしさを感じています。
過去に作ったﾏｸﾛで遅くて話にならないものがあるのですが、
これを利用すれば使い道が出てきそうな気がしてきました。

ﾄﾚﾗﾝｽを導入した事で余計に混乱したのですが、以下は自分用に覚書です。
・分割した空間内に1つしか要素が無い場合、ﾘｽﾄに入れる必要なし
・分割した空間ｻｲｽﾞがﾄﾚﾗﾝｽに満たない場合、正しく処理されない可能性あり

端点一致を探る、組み合わせテスト10

こちらの続きです。
端点一致を探る、組み合わせテスト9 - C#ATIA

昨日の悩みが全て解消するような、良い方法を思い付きました。
時間が無い為、先にｺｰﾄﾞです。

'vba
'ﾓｰﾄﾝ順序を利用した8分木空間分割の為のｸｽﾞｺｰﾄﾞ3
Option Explicit

Private Const CLINER8TREEMANAGER_MAXLEVEL = 7   '有効空間分割最大ﾚﾍﾞﾙ

Sub CATMain()
   Dim min: min = Array(-40#, -40#, -40#)   '最小位置
   Dim Max: Max = Array(40#, 40#, 40#)      '最大位置
   Dim lv&: lv = 2                          '分割ﾚﾍﾞﾙ
   Dim Pnt: Pnt = Array(0, 0, 0)            '調べたい座標値
   Dim Tolerance#: Tolerance = 0.001        '端点一致ﾄﾚﾗﾝｽ
   
   Dim w: w = ArySub(Max, min)
   Dim unit: unit = AryDiv(w, CDbl(sl(1, lv)))
   
   If Not Init(lv) Then Exit Sub
   Dim SpaNo As Variant: SpaNo = GetMortonNum(Pnt, min, unit, Tolerance)
   Dim msg$
   msg = "最小位置 : " + Join(min, ",") + vbNewLine + _
         "最大位置 : " + Join(Max, ",") + vbNewLine + _
         "分割ﾚﾍﾞﾙ : " + CStr(lv) + vbNewLine + _
         "の8分木空間で" + vbNewLine + _
         "座標値  : " + Join(Pnt, ",") + vbNewLine + _
         "線形8分木ｲﾝﾃﾞｯｸｽNo: " + Join(SpaNo, " , ")
   MsgBox msg
End Sub

'*** VBA不足関数 ***

'配列同士の引き算-細かいチェック無し
Private Function ArySub(ByVal a, ByVal b) As Variant
    ArySub = Array(a(0) - b(0), a(1) - b(1), a(2) - b(2))
End Function

'配列と実数の割り算-細かいチェック無し
Private Function AryDiv(ByVal a, ByVal b#) As Variant
    AryDiv = Array(a(0) / b, a(1) / b, a(2) / b)
End Function


'*** Octree ***
'http://marupeke296.com/COL_3D_No15_Octree.html

'線形8分木配列を構築
Private Function Init(ByVal Level&) As Boolean
    Init = False
    
    '設定最高レベル以上の空間は作れない
    If Level >= CLINER8TREEMANAGER_MAXLEVEL Then Exit Function
    
    Init = True
End Function

'ビット分割関数
Private Function BitSeparateFor3D(ByVal n As Byte) As Long
    Dim s As Long: s = n
    s = (s Or sl(s, 8)) And &HF00F
    s = (s Or sl(s, 4)) And &HC30C3
    s = (s Or sl(s, 2)) And &H249249
    BitSeparateFor3D = s
End Function

'8分木モートン順序算出関数
Private Function Get3DMortonNumber(ByVal x As Byte, ByVal y As Byte, ByVal z As Byte) As Long
   Get3DMortonNumber = BitSeparateFor3D(x) Or _
                       sl(BitSeparateFor3D(y), 1) Or _
                       sl(BitSeparateFor3D(z), 2)
End Function

'座標→最小ﾚﾍﾞﾙ空間番号郡取得
''' @param :p-array(Double)-座標値
'''        :Min-array(Double)-空間最小座標値
'''        :unit-array(Double)-空間単位ｻｲｽﾞ
'''        :tol-Double-ﾄﾚﾗﾝｽ
''' @return:array(Long)
Private Function GetMortonNum(ByVal p As Variant, _
                              ByVal min As Variant, _
                              ByVal unit As Variant, _
                              ByVal tol As Double) As Variant
    '最小位置から各軸の比率算出
    Dim ratio#(2), inte&(2), dec#(2) '比率,整数,小数
    Dim tolratio#(2) '各軸のﾄﾚﾗﾝｽの比率
    Dim i&
    For i = 0 To 2
        ratio(i) = (p(i) - min(i)) / unit(i)
        inte(i) = Fix(ratio(i))
        dec(i) = ratio(i) - inte(i)
        tolratio(i) = tol / unit(i)
    Next
    
    'ﾄﾚﾗﾝｽを考慮
    Dim axis(2) As Variant
    Dim spapos As Collection
    For i = 0 To 2
        Set spapos = New Collection
        spapos.add inte(i)
        If dec(i) <= tolratio(i) Then spapos.add inte(i) - 1
        If 1# - dec(i) <= tolratio(i) Then spapos.add inte(i) + 1
        Set axis(i) = spapos
    Next
    
    '空間番号算出
    Dim x, y, z
    Dim SpaNo(): ReDim SpaNo(axis(0).Count * axis(1).Count * axis(2).Count)
    Dim j&, k&, cnt&: cnt = -1
    For Each x In axis(0)
        For Each y In axis(1)
            For Each z In axis(2)
                cnt = cnt + 1
                '本来なら不正値のﾁｪｯｸをすべき
                SpaNo(cnt) = Get3DMortonNumber(x, y, z)
    Next: Next: Next
    ReDim Preserve SpaNo(cnt)
    GetMortonNum = SpaNo
End Function


'*** BitShift ***
'http://www.geocities.co.jp/SiliconValley/4334/unibon/asp/bitshift2.html
' 左シフト
Private Function sl(ByVal x&, ByVal n&) As Long
    If n = 0 Then
        sl = x
    Else
        Dim k: k = CLng(2 ^ (32 - n - 1))
        Dim d: d = x And (k - 1)
        Dim c: c = d * CLng(2 ^ n)
        If x And k Then c = c Or &H80000000
        sl = c
    End If
End Function

細かなﾁｪｯｸ、ｾﾝｽの無い変数・関数名は見逃してください。

今までは点の座標値とﾄﾚﾗﾝｽを考慮し、ﾎﾞﾘｭｰﾑとして線形8分木のｲﾝﾃﾞｯｸｽを
求めていましたが、ﾄﾚﾗﾝｽを考慮して点として最小ﾚﾍﾞﾙのﾓｰﾄﾝ順序の空間番号
を求めるようにしました。(わかり難い表現です･･･)

上記のｺｰﾄﾞを実行した結果はこちら

空間の重心位置(要はど真ん中)の場合、以前ではﾙｰﾄ空間に属する事に
なっていました。　この場合、他に生成された全ての空間と端点一致の
処理が必要でした。　それに対し今回は、最低ﾚﾍﾞﾙの空間の8ヶ所のみの
ﾁｪｯｸで済みます。　この方法であれば
・空間ﾎﾟｲﾝﾀ配列(ｵﾘｼﾞﾅﾙのppCellAry)が、最低ﾚﾍﾞﾙ分のみ
・衝突ﾘｽﾄの巡回が、再帰・非再帰とかではなく単純なﾙｰﾌﾟ1回のみ
・同一空間での要素数が多い際、再度空間分割することが容易
　（極端に離れた要素がある場合への対応が楽）
他にも、ｵﾘｼﾞﾅﾙｺｰﾄﾞではﾁｪｯｸ漏れを起こす可能性がある場合に気が付き
それにも対応したつもりです。

ｵﾘｼﾞﾅﾙｺｰﾄﾞをVBA化する事ばかり考えていたので、計算式の意味を
把握していませんでした。　VBAの商が切り下げじゃないと騒いだ時に
計算式の意味を良く考えれば良かったです。

端点一致を探る、組み合わせテスト9

こちらの続きです。
端点一致を探る、組み合わせテスト8 - C#ATIA

忙しいので進められる状態では無いのですが、考えもまとまらない為
ｺｰﾄﾞ無しです。

最初は端点同士の比較の為、空間を跨ぐ事は無いなぁ　と考えていました。
今回の様に数値を比較する際、ｲｺｰﾙで一致を判断するのはちょっと
まずかったりします。特に実数の場合は。
その為ﾄﾚﾗﾝｽを導入するのですが、この扱いに悩んでます。

組み合わせ数を減らす為空間を分割し、基本的に同一空間と最上位のまでの
親空間との一致を見るだけで済むようになるのですが、極端な表現をすれば
同一ﾚﾍﾞﾙの隣の空間との一致は見ないことになります。
そうなると困るのがﾄﾚﾗﾝｽの扱いでした。

赤いラインは空間の境界ﾗｲﾝです。こんな感じで隣り合った空間に
ﾄﾚﾗﾝｽ以内の距離で端点があると、本来 " 一致する " と判断されるべき
ものが、空間が異なる為に一致を検証されなくなってしまいます。

そこで思いついたのが、端点を重心とし一辺を ﾄﾚﾗﾝｽ X 2 とした
立方体を想定しﾓｰﾄﾝ順序を求めようと考えました。点なのにﾎﾞﾘｭｰﾑで
ﾓｰﾄﾝ順序を求めようとしたのは、これが理由なんです。
本来なら、各空間の境界をﾄﾚﾗﾝｽ分広げてﾓｰﾄﾝ順序を求められると
効率良いだろうと思うのですが、僕の知識ではどの様な計算式にすれば
良いのかわかりませんでした。言葉じゃ表現しにくいのですが･･･。

ﾄﾚﾗﾝｽやﾃﾞｰﾀの状態によって何とも言えませんが、これに引っかかる
要素は恐らく少ない気がするのですが、無視するわけにはいかないので。

さらに問題はここからなんです。　ｵﾘｼﾞﾅﾙの説明での「衝突ｵﾌﾞｼﾞｪｸﾄﾘｽﾄ」と
説明されている部分です。こちらの ④ の説明部分です。
その８ 4分木空間分割を最適化する！

こんなｲﾒｰｼﾞで進めようかとは考えてはいるのですが、

ここで悩むのが、

このチェック後は次の空間に移動する予定なのですが、②の場合は子空間を
持ちません。この場合、スタックにオブジェクトE,Fを登録せずに親の空間に
戻ります。

の部分なんです。
仮に分割ﾚﾍﾞﾙ2とした場合、ﾙｰﾄ空間 - 子空間 - 孫空間　の形になります。
今回は一辺を ﾄﾚﾗﾝｽ X 2 とした立方体なので、要素が極小な為大半は
孫空間に収まるのだろうと想定しています。　その為、生成された孫空間の
上の子空間が生成されない可能性が有り、生成されていない子空間の下の
孫空間はｵﾘｼﾞﾅﾙｻｲﾄの考え方では、巡回されない事になってしまいます。
これは困った･･･。

幸い、ﾙｰﾄ空間・子空間でも自身の下の空間番号がわかる為、最小ﾚﾍﾞﾙ空間
以外の要素を、全て最小ﾚﾍﾞﾙ空間にぶち込んでしまおうかと思ってはいるも
のの、それを行うと一度ﾁｪｯｸした組み合わせを再度ﾁｪｯｸしてしまう事になり、
折角のｱﾙｺﾞﾘｽﾞﾑを効率悪くしてしまいそうな気がして･･･。

しかも、こちらで示した極端に離れた要素がある場合への対処が
上手くできない可能性も。
端点一致を探る、組み合わせテスト4 - C#ATIA

ん～何か良い方法ないかな？