ich habe bei meinem Programm erhebliche Probleme mit der Laufzeit und kenne mich leider noch nicht so gut aus mit MATLAB. Folgendes Problem: Ich habe, wie ihr in meinem Programmcode sehen könnt, mehrere for-Schleifen ineinander verschachtelt. Im Wesentlichen geht es darum, dass die äußerste Schleife, welche sich vermutlich nicht vermeiden lässt, über eine bestimmte Menge von Punkten (ca. 11Mio) einer Punktwolke läuft. Hier werden für jeden Punkt seine k Nachbarn innerhalb eines bestimmten Radius ausgewählt (die entsprechende Index-Suche der Nachbarpunkte wird schon außerhalb der Schleife durchgeführt). Das eigentliche Problem stellen dann aber die folgenden Berechnungen dar, die zwischen jedem Punktpaar innerhalb dieser Nachbarschaft durchgeführt werden müssen. Hierzu verwende ich im Moment zwei weitere for-Schleifen, was den Rechenaufwand bei 11Mio Punkten und jeweils bis zu 30 Nachbarpunkten extrem groß werden lässt.
Ich würde also gerne die beiden inneren for-Schleifen (bzw. wenigstens eine von diesen) ersetzen und weiß nur leider nicht wie...
Hier mein code:
Code:
% aeussere Schleife for i=1:n
... % hier die problematischen Schleifen for b=1:m-1 for a=b+1:m
% Bestimmen von (s)ource und (t)arget if( dot(Pq(a,4:6),Pq(b,1:3)-Pq(a,1:3)) <= ...% na*(pb-pa)
dot(Pq(b,4:6),Pq(a,1:3)-Pq(b,1:3)))% nb*(pa-pb)
ps=Pq(a,1:3); pt=Pq(b,1:3);
ns=Pq(a,4:6); nt=Pq(b,4:6);
else
ps=Pq(b,1:3); pt=Pq(a,1:3);
ns=Pq(b,4:6); nt=Pq(a,4:6);
end % Darboux-Rahmen (uvw)
u=ns; v=cross(pt-ps,u)/norm(pt-ps); w=cross(u,v);
% Berechnen der Histogramm-Features f0, f1 und f2
f0 = dot(v,nt);
f1 = dot(u,pt-ps)/norm(pt-ps);
f2 = atan2(dot(w,nt),dot(u,nt)) /pi;
idx = step(d, [f0 f1 f2]);
PFHq(idx) = PFHq(idx) + 1/anzPP;
end end ... end
Wäre schön, wenn ihr Vorschläge habt, wie ich die Berechnung optimieren kann! Falls ich was Wichtiges vergessen habe, oder ihr zusätzliche Informationen braucht, dann gebt Bescheid!
Schon vielen Dank für eure Mühe!
Ja, das habe ich eigentlich schon berücksichtigt. Also alle Matrizen, in denen ich Werte abspeichere, habe ich vorher mit der "zeros"-Funktion erzeugt.
Hat noch keiner eine Idee? Ist die Berechnung vielleicht gar nicht anders möglich? Ich habe mir selbst natürlich auch schon einige Gedanken gemacht und ich habe überlegt die Indizes der Punktpaare, welche durch die for-Schleifen erzeugt werden einfach in Vektor-/Matrizenform zu verwenden und somit die for-Schleifen zu umgehen. Ich weiß allerdings nicht, wie ich das dann auf die Berechnungen übertragen soll, die zwischen jedem einzelnen Punktpaar durchgeführt werden müssen.
Gut zu wissen, dann lege ich gerne mal nach!
Hier also die ausführlichere (hoffentlich verständliche) Erklärung der Berechungen und der ausführlichere code:
Zunächst sind die Eingangsdaten für das Programm in der Matrix "P" gespeichert. Dabei handelt es sich um eine [~11Mio. x 3]-Matrix mit den dreidimensionalen Koordinaten von Scanpunkten. Dann werden die Flächennormalen zu allen Punkten berechnet und zusammen mit den Koordinaten in einer [~11Mio. x 6]-Matrix gespeichert: PN = [X Y Z Nx Ny Nz]
(den code hierzu habe ich mal weggelassen!)
Anschließend soll für jeden Punkt ein sogenanntes "Point Feature Histograms (PFH)" berechnet werden (im folgenden code zu sehen).
Das PFH eines Abfragepunktes pq ist schließlich ein [1 x 27]-Vektor, mit positiven Zahlen, die aufsummiert "1" ergeben.
Dazu werden zuerst die "k" Nachbarpunkte jeden Punktes innerhalb eines Radius gesucht, bzw. deren Indizes.
Bei der PFH-Berechnung läuft dann die äußerste Schleife über alle 11Mio. Punkte. Für jeden Punkt werden dann zunächst seine Nachbarpunkte ausgewählt und in der [k x 6]-Matrix "Pq" gespeichert. Sie enthält also untereinander den aktuellen Punkte und die zugehörigen Nachbarpunkte.
Aus dieser Nachbarschaft werden dann zwischen jedem Paar von Punkten Berechnungen durchgeführt. Um jedes Punktpaar auszuwählen dienen die beiden (problematischen) verschachtelten for-Schleifen.
Innerhalb dieser Schleifen wird zuerst bestimmt welcher der beiden aktuellen Punkte der source- und welcher der target-point ist. Dann wird ein Darboux-Rahmen mit den [1 x 3]-Vektoren "u", "v" und "w" aufgestellt und daraus die Features "f0", "f1" und "f2" berechnet, die die Oberfläche um den Abfragepunkt beschreiben. Bei diesen Features handelt es sich um skalare Größen, deren Werte zwischen -1 und 1 liegen.
Dann gibt die (eigene) Funktion "step" einen Index {1,2,...27} "idx" zurück, welcher sich aus den zuvor berechneten Features ergibt. Dieser Index gibt an, an welcher Stelle des PFHs das aktuelle Punktpaar liegt und speichert den entsprechenden prozentualen Anteil an der Stelle "idx" in "PFHq".
Die umschließende "if-Abfrage" dient der Überprüfung, ob genügend Nachbarpunkte für die Berechnung vorhanden sind und ob die Flächennormale des Abfragepunktes vorhanden ist.
Code:
%% Finden der k-Nachbarschaft (PFH-Berechnung) % Radius für Nachbarschaftssuche
r_PFH = 2.5; %[m] % minimale / maximale(k) #Nachbarn
min_k_PFH = 5;
k_PFH = 21;
% Radius-Suche
IDX_PFH = rangesearch(PN(:,1:3),PN(:,1:3),r_PFH);
%% Point Feature Histogram (PFH) % #Wert-Intervalle pro Feature
d = 3; % d^3=#Bins im Histogramm % maximale #Punktpaare pro Abfragepunkt pq
anzPP_max = k_PFH*(k_PFH+1)/2;
% Initialisieren der PFH-Matrix [n x #Bins] mit n...Anzahl Punkte (~11Mio.)
PFH = zeros(n,d^3);
% Start der for-Schleife zur Berechnung der PFHs for i=1:n
% Ueberpruefung ob Voraussetzung fuer PFH-Berechnung erfuellt ist if(length(IDX_PFH{i}) >= min_k_PFH+1 &&... % genuegend Nachbarn vorh. sum(PN(i,4:6)) ~= 0)% Normale existiert % Abfragepunkt pq mit k-Nachbarschaft if(length(IDX_PFH{i}) <= k_PFH+1)% keine Ueberschreitung von k
Pq = PN(IDX_PFH{i},:);
else
Pq = PN(IDX_PFH{i}(1:k_PFH+1),:);
end % #Punkte
m = length(Pq);
% #Punktpaare
anzPP = (m-1)*m/2;
% Initialisieren des PFHs für pq
PFHq = zeros(1,d^3);
for b=1:m-1 for a=b+1:m
% Bestimmen von (s)ource und (t)arget if( dot(Pq(a,4:6),Pq(b,1:3)-Pq(a,1:3)) <= ...% na*(pb-pa)
dot(Pq(b,4:6),Pq(a,1:3)-Pq(b,1:3)))% nb*(pa-pb)
ps=Pq(a,1:3); pt=Pq(b,1:3);
ns=Pq(a,4:6); nt=Pq(b,4:6);
else
ps=Pq(b,1:3); pt=Pq(a,1:3);
ns=Pq(b,4:6); nt=Pq(a,4:6);
end % Darboux-Rahmen (uvw)
u=ns; v=cross(pt-ps,u)/norm(pt-ps); w=cross(u,v);
% Berechnen der Histogramm-Features f0, f1 und f2
f0 = dot(v,nt);
f1 = dot(u,pt-ps)/norm(pt-ps);
f2 = atan2(dot(w,nt),dot(u,nt)) /pi; % -pi...pi => -1...1 % Zuweisen des Histogram-Index ueber berechnete Features
idx = step(d, [f0 f1 f2]);
% Aktualisieren des PFHs
PFHq(idx) = PFHq(idx) + 1/anzPP;
end end % Abspeichern des PFHs
PFH(i,:) = PFHq;
else% Nullsetzen des i-ten PFHs, falls nicht genuegend Nachbarn vorh.
PFH(i,:) = zeros(1,d^3);
end end
Ich habe mir selbst natürlich auch schon einige Gedanken gemacht und ich habe überlegt die Indizes der Punktpaare, welche durch die for-Schleifen erzeugt werden einfach in Vektor-/Matrizenform zu verwenden und somit die for-Schleifen zu umgehen. Ich weiß allerdings nicht, wie ich das dann auf die Berechnungen übertragen soll, die zwischen jedem einzelnen Punktpaar durchgeführt werden müssen.
Zu meiner Idee mit den Indizes habe ich mal ein kleines Programm geschrieben. Hier werden einfach ein paar zufällige Punkte erzeugt und dann für jedes Punktpaar jeweils der Mittelwert berechnet.
Hier der lauffähige code:
Code:
%% for-Schleifen - TEST % Dieses Testprogramm erzeugt zunaechst eine Matrix mit zufaelligen % Testpunkten. Ziel ist fuer jedes Paar von Testpunkten den zugehoerigen % Mittelwert zu berechnen. % Drei verschiedene Varianten werden verwendet, die auf dasselbe Ergebnis % fuehren: Verwendung von zwei, einer bzw. keiner for-Schleife/n
%% Testdaten und Initialisierung % #Punkte
n = 5;
% Erzeugen von zufaelligen Punkten (Werte von 0 bis 4)
P = floor(5.*rand(n,3)); % [n x 3] % #Punktpaare
anzPP = (n-1)*n/2;
% Initialisieren der Ergebnismatrizen
M1 = zeros(anzPP,3); M2 = zeros(anzPP,3);
% Initialisieren der Indexmatrizen
I = zeros(n-1,n-1); J=I;
%% Variante 1: Zwei for-Schleifen % Laufvariable fuer Zuordnung
z=1;
for i=1:n-1 for j=i+1:n
% Indexmatrizen fuer weitere Berechnungen (M2 und M3)
I(i,j-1)=i;
J(i,j-1)=j;
% Berechnung des Mittelwertes und Speichern in Ergebnismatrix M1
M1(z,:) = (P(i,:) + P(j,:))/2;
z=z+1;
end end
%% Variante 2: Eine for-Schleife % Laufvariablen fuer Zuordnung
a=1; b=n-1; c=n-2;
for d=1:n-1 % Berechnung der Mittelwerte und Speichern in Ergebnismatrix M2
M2(a:b,:) = (P(I(d,d:end),:) + P(J(d,d:end),:))/2;
a=b+1; b=b+c; c=c-1;
end
%% Variante 3: Keine for-Schleife % Umformen der Indexmatrizen
I2 = reshape(I',(n-1)^2,1); I2 = I2(I2~=0);
J2 = reshape(J',(n-1)^2,1); J2 = J2(J2~=0);
% Berechnung der Mittelwerte und Speichern in Ergebnismatrix M2
M3 = (P(I2,:) + P(J2,:))/2;
Hier kann ich also problemlos die zwei verschachtelten for-Schleifen umgehen, indem ich die Berechnung über die Indizes der Punktpaare durchführe.
Aber meine Frage ist nun, ob (und falls ja, wie??) ich diese Idee auf mein eigentlich Problem anwenden kann?
Einstellungen und Berechtigungen
Du kannst Beiträge in dieses Forum schreiben. Du kannst auf Beiträge in diesem Forum antworten. Du kannst deine Beiträge in diesem Forum nicht bearbeiten. Du kannst deine Beiträge in diesem Forum nicht löschen. Du kannst an Umfragen in diesem Forum nicht mitmachen. Du kannst Dateien in diesem Forum posten Du kannst Dateien in diesem Forum herunterladen
MATLAB, Simulink, Stateflow, Handle Graphics, Real-Time Workshop, SimBiology, SimHydraulics, SimEvents, and xPC TargetBox are registered trademarks and The MathWorks, the L-shaped membrane logo, and Embedded MATLAB are trademarks of The MathWorks, Inc.
RSS