Abt. Fast Bounding Sphere
=================
Mit einem in früheren Beitragskolumnen veröffentlichten Programm "Fliegenschwärme fangen" haben wir die 'Konvexe Hülle' einer Punktwolke halbwegs flott dargestellt. Wollten wir so einen Schwarm aber tatsächlich fangen, so könnten wir nur ein Schmetterlingsnetz benutzen. Schmetterlingsnetze sind üblicherweise auf einen runden Metallrahmen gespannt. Nur: Wie groß müsste der bei einer bestimmten Punktwolke mindestens sein und wohin würden wir ungefähr zielen?
Ein Algorithmus, der hier zu genau rechnet, verpasst mit Sicherheit einige Fliegen, weil die sich ja ziemlich flott dahinbewegen. Wenn eine ausreichend gute und schnelle Berechnungsmethode gesucht ist, bietet sich der Algorithmus von [Ritter 1990] an, der nachstehend in XProfan übersetzt wurde.
Gruss
P.S.: Die Originalvorlage benutzt Vektorarithmetik und kann n-dimensionale Fliegenschwärme fangen, für Demozwecke habe ich das in XProfan aber auf 2D-Punktwolken begrenzt. Eine Realisation für mehr Dimensionen sollte aber eigentlich keine Hexerei sein.
WindowTitle "Ritter's Fast Bounding-Ball Algorithm "+\
"(Errechnet eine ziemlich kleine Punktwolken-umfassende Kugel, hier in 2D,"+\
" im Schnitt 5% bis 20% größer als Minimalsphäre, aber sehr schnell!"
' <TransCoded to XProfan 11, (CL) CopyLeft 2014-07 by P. Specht@gmx.at>
' <Source: http:'geomalgorithms.com/a08-_containers.html#fastBall%28%29>
' Theory: http://en.wikipedia.org/wiki/Bounding_sphere
' =====================================================================
' Copyright 2001 softSurfer, 2012 Dan Sunday
' This code may be freely used and modified for any purpose providing that this
' copyright notice is included with it. SoftSurfer makes no warranty for this code,
' and cannot be held liable for any real or imagined damage resulting from its use.
' Users of this code must verify correctness for their application.
' <The Program is> based on the algorithm given by [Jack Ritter, 1990].
' ====================================================================================
' <ASSUMPTIONS ON AVAILABLE VECTOR MATH OVERRIDDEN BY SIMPLE XPROFAN MATH>
' Result: Ball (2D: Circle) B with Center and Radius <Remarks <..> by P.Specht>
'=====================================================================================
' TESTSYSTEM mit n& Punkten:
var n&=12:dec n&
'=====================================================================================
DECLARE BCenterX!,BCentery!,BRadius! ' Übergabevariablen, global definiert
Declare Px![n&],Py![n&],i&
WindowStyle 24:Font 2:Randomize:Window 0,0-%maxx,%maxy-40
var xh&=width(%hwnd)/2:var yh&=height(%hwnd)/2
REPEAT
Px![]=xh&*0.75+rnd(xh&/2):Py![]=yh&*1.5-rnd(xh&/2)
CLS:usepen 0,5,255:whileloop 0,n&:i&=&Loop:line Px![i&],Py![i&] - Px![i&]+1,Py![i&]:endwhile
FastBall Px![],Py![],n&
Usepen 0,2,rgb(0,0,255):usebrush 0,rgb(0,255,0)
Ellipse BCenterX!-BRadius!,BCenterY!+BRadius! - BCenterX!+BRadius!,BCenterY!-BRadius!
waitinput 10000
until %key = 27
END
proc FastBall
parameters Px![],Py![],n&
' (based on the algorithm given by [Jack Ritter, 1990])
declare Cx!,Cy!,rad!,rad2! ' Center of ball, radius and radius squared
declare xmin!, xmax!, ymin!, ymax! ' bounding box extremes
declare Pxmin&,Pxmax&,Pymin&,Pymax& ' index of P[] at box extreme
declare i&,dPxx!,dPxy!,dPyx!,dPyy!,dx2!,dy2!,dist!,dist2!,dPx!,dPy!
' Find a large diameter to start with
xmin!=10^50:xmax!=-1*10^50:ymin!=10^50:ymax!=-1*10^50:Pxmin&=0:Pxmax&=0:Pymin&=0:Pymax&=0
whileloop 0,n&:i&=&Loop
if Px![i&]<xmin!: xmin!=Px![i&]:Pxmin&=i&:endif
if Px![i&]>xmax!: xmax!=Px![i&]:Pxmax&=i&:endif
if Py![i&]<ymin!: ymin!=Py![i&]:Pymin&=i&:endif
if Py![i&]>ymax!: ymax!=Py![i&]:Pymax&=i&:endif
endwhile
:::::::Usepen 0,1,rgb(0,0,255):usebrush 0,rgb(0,255,0)
:::::::rectangle Px![Pxmin&],Py![Pymin&] - Px![Pxmax&],Py![Pymax&]
cx!=(Px![Pxmin&]+Px![Pxmax&])/2
cy!=(Py![Pymin&]+Py![Pymax&])/2
rad2!=sqr(cx!-Px![Pxmin&])
case rad2!<sqr(cy!-Py![Pymin&]):rad2!=sqr(cy!-Py![Pymin&])
rad!=sqrt(rad2!)
:::::::usepen 0,9,0:line cx!,cy! - cx!+1,cy!
' Now check that all points P[i] are inside the Ball
' and if not, expand the ball just enough to include them
' Vector dP; float dist, dist2;
whileloop 0,n&:i&=&Loop
dPx!=Px![i&]-Cx!
dPy!=Py![i&]-Cy!
dist2!=sqr(dPx!)+sqr(dPy!) '= norm2(dP);
case dist2!<=rad2!:continue ' P[i] is inside the ball already
' P[i] not in ball, so expand ball to include it:
dist!=sqrt(dist2!)
rad!=(rad!+dist!)/2 ' enlarge radius just enough
rad2!=sqr(rad!)
Cx!=Cx!+((dist!-rad!)/dist!)*dPx! ' shift Center toward P[i]
Cy!=Cy!+((dist!-rad!)/dist!)*dPy!
Endwhile
Bcenterx!=Cx!:Bcentery!=Cy!:BRadius!=rad!
Endproc
Alles anzeigen