C#ATIA

↑タイトル詐欺 主にFusion360API 偶にCATIA V5 VBA(絶賛ネタ切れ中)

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

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

一応、取り組んでますが、全くCATIAとは無関係な状態です。

'vba
'モートン順序を利用した8分木空間分割の為のクズコード
Option Explicit

Sub CATMain()
   Dim min: min = Array(-12#, -34#, -56#)   '最小位置
   Dim max: max = Array(98#, 76#, 54#)      '最大位置
   Dim lv&: lv = 2                          '分割レベル
   Dim pnt: pnt = Array(12.3, 45.6, 7.89)   '調べたい座標値
   
   Dim w: w = ArySub(max, min)
   Dim unit: unit = AryDiv(w, CDbl(sl(1, lv)))
   Dim morton&: morton = GetPointElem(pnt, min, unit)
   
   Dim msg$
   msg = "最小位置 : " + Join(min, ",") + vbNewLine + _
         "最大位置 : " + Join(max, ",") + vbNewLine + _
         "分割レベル : " + CStr(lv) + vbNewLine + _
         "の8分木空間で" + vbNewLine + _
         "座標値   : " + Join(pnt, ",") + vbNewLine + _
         "のモートン順序は [ " + CStr(morton) + " ] です"
   MsgBox msg
End Sub

'配列同士の引き算-細かいチェック無し
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

'ビット分割関数
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

'座標→線形8分木要素番号変換関数
''' @param :p-array(Double)-変更するかも
'''        :rgnMin-array(Double)-変更するかも
'''        :unit-array(Double)-変更するかも
''' @return:Long
Private Function GetPointElem(ByVal p As Variant, _
                              ByVal rgnMin As Variant, _
                              ByVal unit As Variant) As Long
    GetPointElem = Get3DMortonNumber( _
            CByte((p(0) - rgnMin(0)) / unit(0)), _
            CByte((p(1) - rgnMin(1)) / unit(1)), _
            CByte((p(2) - rgnMin(2)) / unit(2)))
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

' 右シフト(算術(>>)ではなく論理(>>>)シフトに相当)
Private Function sr(ByVal x&, ByVal n&) As Long
    If n = 0 Then
        sr = x
    Else
        Dim y: y = x And &H7FFFFFFF
        Dim z
        If n = 32 - 1 Then
            z = 0
        Else
            z = y \ CLng(2 ^ n)
        End If
        If y <> x Then z = z Or CLng(2 ^ (32 - n - 1))
        sr = z
    End If
End Function

' オーバフローを無視して 32 ビットの加算をおこなう
Private Function add(ByVal a&, ByVal b&) As Long
    If a >= 0 And b <= 0 Then
        add = a + b
    ElseIf a <= 0 And b >= 0 Then
        add = a + b
    Else
        Dim x: x = a And &H3FFFFFFF
        Dim y: y = b And &H3FFFFFFF
        Dim z: z = x + y
        Dim f: f = 0
        
        If z And &H40000000 Then f = f + 1
        z = z And &H3FFFFFFF
        If a And &H40000000 Then f = f + 1
        If a And &H80000000 Then f = f + 2
        If b And &H40000000 Then f = f + 1
        If b And &H80000000 Then f = f + 2
        If f And 1 Then z = z Or &H40000000
        If f And 2 Then z = z Or &H80000000
        add = z
    End If
End Function

'オーバフローを無視して 32 ビットの加算をおこなう
Private Function addCur(ByVal a&, ByVal b&) As Long
    Dim c@: c = CCur(a) + CCur(b)
    If c > &H7FFFFFFF Then
        c = c - CCur(2 ^ 32)
    ElseIf c < &H80000000 Then
        c = c + CCur(2 ^ 32)
    End If
    addCur = CLng(c)
End Function

ビットシフトについては、こちらのサイトを参考に型定義した形に
修正しました。(ありがとうございます)
VBScript で高速にビットシフト演算する (unibon)
左シフトのsl関数しか利用しておりませんが・・・。

座標値からモートン順序を算出する、GetPointElem関数はちょっと考えがあって
同一ではありません。

上記コードの実行結果はこちら
f:id:kandennti:20161213173831p:plain
正直に書くと、正しいのか確認しておりません。
(間違いを指摘していただけると助かります)