Converting a short piece of C/C++ code

[IdeasVacuum]

I found that not paying attention to the plane normal direction had an adverse result, as indeed did the relative success of creating planes - in some cases, planes were being created correctly but the point sets they were based on required better filtering. So, with "better" planes comes better intersection lines.

However, I still don't have a good 3D Point at Line-Line intersection method.

[alter Mann]

Try this:

Code: Select all

Procedure DistancePointLine ( Array adQ.d(1), Array adP.d(1), Array adV.d(1), *dU.double, *dA.double)
  Protected Dim adU.d(2)
  Protected Dim adN.d(2)
  Protected dL
  
  dL = adV(0)*adV(0) + adV(1)*adV(1) + adV(2)*adV(2)
  *dA\d = 1.0e20
  *dU\d = 0.0
  
  adU(0) = adQ(0) - adP(0)
  adU(1) = adQ(1) - adP(1)
  adU(2) = adQ(2) - adP(2)
  
  If dL > 0.0
    adN(0) = adU(1)*adV(2) - adU(2)*adV(1)
    adN(1) = adU(2)*adV(0) - adU(0)*adV(2)
    adN(2) = adU(0)*adV(1) - adU(1)*adV(0)
    *dA\d = (adN(0)*adN(0) + adN(1)*adN(1) + adN(2)*adN(2)) / dL
    *dU\d = (adU(0)*adV(0) + adU(1)*adV(1) + adU(2)*adV(2)) / dL
  Else
    *dA\d = adU(0)*adU(0) + adU(1)*adU(1) + adU(2)*adU(2)
  EndIf
  
  *dA\d = Sqr(*dA\d)
EndProcedure

Procedure.d VecSpat (Array adU.d(1), Array adV.d(1), Array adW.d(1))
  ProcedureReturn adU(0)*adV(1)*adW(2)+adV(0)*adW(1)*adU(2)+adW(0)*adU(1)*adV(2) - adU(0)*adW(1)*adV(2)-adV(0)*adU(1)*adW(2)-adW(0)*adV(1)*adU(2)
EndProcedure

Procedure DistanceLineLine (Array adP.d(1), Array adU.d(1), Array adQ.d(1), Array adV.d(1), *dU.double, *dV.double, *dA.double)
;adP - point of first line
;adV - direction of first line
;adQ - point of second line
;adV - direction of second line
;results : *dU - parameter of shortest distance for first line - to calculate point use   adP + *dU\d * adU
;             *dV - parameter of shortest distance for second line - to calculate point use   adQ + *dV\d * adV
;             *dA - shortest distance between the two lines   - *dA\d=0 -> intersection
  Protected Dim adW.d(2)
  Protected Dim adX.d(2)
  Protected dL.d
  adW(0) = adU(1)*adV(2) - adU(2)*adV(1)
  adW(1) = adU(2)*adV(0) - adU(0)*adV(2)
  adW(2) = adU(0)*adV(1) - adU(1)*adV(0)
  dL = Sqr(adW(0)*adW(0)+adW(1)*adW(1)+adW(2)*adW(2))
  If dL > 1.0e-10 ; lines are not colinear
    dL = 1.0 / dL
    adW(0) * dL
    adW(1) * dL
    adW(2) * dL
    dL = 1.0/VecSpat(adU(),adV(),adW())
    adX(0) = adQ(0)-adP(0)
    adX(1) = adQ(1)-adP(1)
    adX(2) = adQ(2)-adP(2)
    *dU\d = VecSpat(adX(),adV(),adW()) * dL
    *dV\d = VecSpat(adU(),adX(),adW()) * dL
    *dA\d = Abs(VecSpat(adU(),adV(),adX())*dL)
  ElseIf adU(0)*adU(0)+adU(1)*adU(1)+adU(2)*adU(2) > 1.0e-20
    DistancePointLine(adQ(),adP(),adU(),*dU,*dA)
    If adV(0)*adV(0)+adV(1)*adV(1)+adV(2)*adV(2) > 1.0e-20
      DistancePointLine(adP(),adQ(),adV(),*dV,*dA)
    Else
      *dV\d = 0.0
    EndIf
  ElseIf adV(0)*adV(0)+adV(1)*adV(1)+adV(2)*adV(2) > 1.0e-20
    *dU\d = 0.0
    DistancePointLine(adP(),adQ(),adV(),*dV,*dA)
  Else
    *dU\d = 0.0
    *dV\d = 0.0
    adW(0) = adQ(0)-adP(0)
    adW(1) = adQ(1)-adP(1)
    adW(2) = adQ(2)-adP(2)
    *dA\d = Sqr(adW(0)*adW(0)+adW(1)*adW(1)+adW(2)*adW(2))
  EndIf
EndProcedure

[IdeasVacuum]

Thanks alter Mann - I'm going to try your code with stronger coffee

[alter Mann]

f you are sure that the lines are not parallel, you can shorten the code:

Code: Select all

Procedure.d VecSpat (Array adU.d(1), Array adV.d(1), Array adW.d(1))
  ProcedureReturn adU(0)*adV(1)*adW(2)+adV(0)*adW(1)*adU(2)+adW(0)*adU(1)*adV(2) - adU(0)*adW(1)*adV(2)-adV(0)*adU(1)*adW(2)-adW(0)*adV(1)*adU(2)
EndProcedure

Procedure DistanceLineLine (Array adP.d(1), Array adU.d(1), Array adQ.d(1), Array adV.d(1), *dU.double, *dV.double, *dA.double)
;adP - point of first line
;adV - direction of first line
;adQ - point of second line
;adV - direction of second line
;results : *dU - parameter of shortest distance for first line - to calculate point use   adP + *dU\d * adU
;             *dV - parameter of shortest distance for second line - to calculate point use   adQ + *dV\d * adV
;             *dA - shortest distance between the two lines   - *dA\d=0 -> intersection
  Protected Dim adW.d(2)
  Protected Dim adX.d(2)
  Protected dL.d
  adW(0) = adU(1)*adV(2) - adU(2)*adV(1)
  adW(1) = adU(2)*adV(0) - adU(0)*adV(2)
  adW(2) = adU(0)*adV(1) - adU(1)*adV(0)
  dL = Sqr(adW(0)*adW(0)+adW(1)*adW(1)+adW(2)*adW(2))
  If dL > 1.0e-10 ; lines are not colinear
    dL = 1.0 / dL
    adW(0) * dL
    adW(1) * dL
    adW(2) * dL
    dL = 1.0/VecSpat(adU(),adV(),adW())
    adX(0) = adQ(0)-adP(0)
    adX(1) = adQ(1)-adP(1)
    adX(2) = adQ(2)-adP(2)
    *dU\d = VecSpat(adX(),adV(),adW()) * dL
    *dV\d = VecSpat(adU(),adX(),adW()) * dL
    *dA\d = Abs(VecSpat(adU(),adV(),adX())*dL)
    ProcedureReturn #True ; nearest distance was found
  EndIf
  ProcedureReturn #False
EndProcedure

[IdeasVacuum]

.... yes, for further-down-the-process reasons the lines get selected in anti-clockwise order, so two lines parallel would not occur.