2013-10-29 18 views
9

Sto cercando di trovare l'equazione parametrica della traiettoria di un punto di partenza su diversi punti sulla superficie di una sfera unitaria tale che:traiettoria ottimale sulla superficie di una sfera

  1. ogni salto è piccolo (pi/4 < d < pi/2) e in un intervallo ristretto, ad es. [1.33, 1.34]
  2. le visite punto più regioni della sfera modo più rapido e uniforme possibile
  3. il punto viaggia lungo "vettori di direzione" più diverso possibile

Questo è quello che ho provato

N = 3600; % number of points 
t = (1:N) * pi/180; % parameter 
theta_sph = sqrt(2) * t * pi; % first angle 
phi_sph = sqrt(3) * t * pi; % second angle 
rho_sph = 1; % radius 
% Coordinates of a point on the surface of a sphere 
x_sph = rho_sph * sin(phi_sph) .* cos(theta_sph); 
y_sph = rho_sph * sin(phi_sph) .* sin(theta_sph); 
z_sph = rho_sph * cos(phi_sph); 

% Check length of jumps (it is intended that this is valid only for small jumps!!!) 
aa = [x_sph(1:(N-1)); y_sph(1:(N-1)); z_sph(1:(N-1))]; 
bb = [x_sph(2:N); y_sph(2:N); z_sph(2:N)]; 
cc = cross(aa, bb); 
d = rho_sph * atan2(arrayfun(@(n) norm(cc(:, n)), 1:size(cc,2)), dot(aa, bb)); 
figure 
plot(d, '.') 
figure 
plot(diff(d), '.') 

% Check trajectory on the surface of the sphere 
figure 
hh = 1; 
h_plot3 = plot3(x_sph(hh), y_sph(hh), z_sph(hh), '-'); 
hold on 
axis square 
% axis off 
set(gca, 'XLim', [-1 1]) 
set(gca, 'YLim', [-1 1]) 
set(gca, 'ZLim', [-1 1]) 
for hh = 1:N 
    h_point3 = plot3(x_sph(hh), y_sph(hh), z_sph(hh), ... 
     'o', 'MarkerFaceColor', 'r', 'MarkerEdgeColor', 'r'); 
    drawnow 
    delete(h_point3) 
    set(h_plot3, 'XData', x_sph(1:hh)) 
    set(h_plot3, 'YData', y_sph(1:hh)) 
    set(h_plot3, 'ZData', z_sph(1:hh)) 
end 

EDIT -> qualcuno può trovare una traiettoria più regolare, forse che copre la sfera più velocemente (cioè con il minor numero di salti) e più uniforme? Traiettoria regolare, nel senso che dovrebbe cambiare direzione senza intoppi, non bruscamente. La bellezza estetica è un bonus. I punti dovrebbero essere distribuiti sulla superficie della sfera nel modo più uniforme possibile.

+1

Potete chiarire ulteriormente i criteri? In particolare i punti 2 e 3, e anche cosa intendi per "normale" e "uniforme" alla fine? – aganders3

+0

Uniforme significa che i punti dovrebbero essere sparsi su tutta la superficie della sfera, i più lontani tra loro (questo può essere misurato bene). In modo improprio, ovviamente, proprio come un'intuizione generica: per traiettoria normale intendevo una traiettoria che non cambia troppo, cioè dovrebbe cambiare direzione senza intoppi, non bruscamente. E se è esteticamente gradevole, questo è un bonus :-) –

+1

Puoi chiarire il punto 3?Se lo ignoro, mi aspetto che una spirale lungo la sfera (che è il percorso che stai utilizzando) sia ottimale. Basta regolare il numero di giri completi che l'angolo polare attraversa per controllare la dimensione e la copertura del salto. Se controlli la velocità polare puoi "barcollare" i tuoi punti in modo che non cadano tutti sulle stesse linee di latitudine poiché sembra che tu sia più preoccupato per la copertura del "punto" rispetto alla traiettoria continua. – kalhartt

risposta

2

Modifica dell'inizio del codice per introdurre una rotazione della sfera sottostante. Questo dà una traiettoria che non ritorna più frequentemente ai poli. Potrebbe richiedere un po 'di messa a punto delle velocità di rotazione per sembrare "bello" (e probabilmente sembra migliore quando ruota solo attorno ad un asse, non tutti i 3). rot_angle1 è la rotazione attorno all'asse x e rot_angle2 e rot_angle3 sono la rotazione attorno agli assi y e z. Forse questo ti dà almeno un'idea!

N = 3600; % number of points 
t = (1:N) * pi/180; % parameter 
theta_sph = sqrt(2) * t * pi; % first angle 
phi_sph = sqrt(3) * t * pi; % second angle 
rho_sph = 1; % radius 
rot_angle1 = sqrt(2) * t * pi; 
rot_angle2 = sqrt(2.5) * t * pi; 
rot_angle3 = sqrt(3) * t * pi; 
% Coordinates of a point on the surface of a sphere 
x_sph0 = rho_sph * sin(phi_sph) .* cos(theta_sph); 
y_sph0 = rho_sph * sin(phi_sph) .* sin(theta_sph); 
z_sph0 = rho_sph * cos(phi_sph); 

x_sph1 = x_sph0; 
y_sph1 = y_sph0.*cos(rot_angle1)-z_sph0.*sin(rot_angle1); 
z_sph1 = y_sph0.*sin(rot_angle1)+z_sph0.*cos(rot_angle1); 

x_sph2 = x_sph1.*cos(rot_angle2)+z_sph1.*sin(rot_angle2); 
y_sph2 = y_sph1; 
z_sph2 = -x_sph1.*sin(rot_angle2)+z_sph1.*cos(rot_angle2); 

x_sph = x_sph2.*cos(rot_angle3)-y_sph2.*sin(rot_angle3); 
y_sph = x_sph2.*sin(rot_angle3)+y_sph2.*cos(rot_angle3); 
z_sph = z_sph2; 
0

modifico qui perché è un codice lungo. Dopo i suggerimenti di David e Kalhartt, ho provato questo:

N = 3600; % number of points 
t = (1:N) * pi/180; % parameter 
% theta_sph much faster than phi_sph to avoid overly visiting the poles 
theta_sph = sqrt(20.01) * t * pi; % first angle 
phi_sph = sqrt(.02) * t * pi; % second angle 
rho_sph = 1; % radius 
% Coordinates of a point on the surface of a sphere 
x_sph0 = rho_sph * sin(phi_sph) .* cos(theta_sph); 
y_sph0 = rho_sph * sin(phi_sph) .* sin(theta_sph); 
z_sph0 = rho_sph * cos(phi_sph); 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
% Use David's hint to rotate axes (but only the first now): 
rot_angle1 = t * pi/10; 
rot_angle2 = 0; 
rot_angle3 = 0; 

x_sph1 = x_sph0; 
y_sph1 = y_sph0.*cos(rot_angle1)-z_sph0.*sin(rot_angle1); 
z_sph1 = y_sph0.*sin(rot_angle1)+z_sph0.*cos(rot_angle1); 

x_sph2 = x_sph1.*cos(rot_angle2)+z_sph1.*sin(rot_angle2); 
y_sph2 = y_sph1; 
z_sph2 = -x_sph1.*sin(rot_angle2)+z_sph1.*cos(rot_angle2); 

x_sph = x_sph2.*cos(rot_angle3)-y_sph2.*sin(rot_angle3); 
y_sph = x_sph2.*sin(rot_angle3)+y_sph2.*cos(rot_angle3); 
z_sph = z_sph2; 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 

% Check length of jumps 
aa = [x_sph(1:(N-1)); y_sph(1:(N-1)); z_sph(1:(N-1))]; 
bb = [x_sph(2:N); y_sph(2:N); z_sph(2:N)]; 
cc = cross(aa, bb); 
d = rho_sph * atan2(arrayfun(@(n) norm(cc(:, n)), 1:size(cc,2)), dot(aa, bb)); 
figure 
plot(d, '.') 

% Check trajectory on the surface of the sphere 
figure 
hh = 1; 
h_plot3 = plot3(x_sph(hh), y_sph(hh), z_sph(hh), '-'); 
hold on 
axis square 
% axis off 
set(gca, 'XLim', [-1 1]) 
set(gca, 'YLim', [-1 1]) 
set(gca, 'ZLim', [-1 1]) 
for hh = 1:N 
    h_point3 = plot3(x_sph(hh), y_sph(hh), z_sph(hh), ... 
     'o', 'MarkerFaceColor', 'r', 'MarkerEdgeColor', 'r'); 
    drawnow 
    delete(h_point3) 
    set(h_plot3, 'XData', x_sph(1:hh)) 
    set(h_plot3, 'YData', y_sph(1:hh)) 
    set(h_plot3, 'ZData', z_sph(1:hh)) 
end 

Penso che sia molto meglio di prima! Due cose che ho trovato importanti: 1. theta_sph deve essere molto più veloce di phi_sph per evitare di visitare due poli opposti troppo spesso; 2. se theta_sph va più veloce di phi_sph, allora devi ruotare lentamente su rot_angle1 o rot_angle2 per ottenere una traiettoria che non sia troppo disordinata. Sono ancora aperto a qualsiasi altro suggerimento per migliorare il risultato.

1

Non ho una copia di matlab a portata di mano, ma posterò le modifiche che apporterei alla tua curva.

Giusto per essere chiari poiché ci sono n-finity + 1 diverse definizioni per gli angoli sferici. Userò quanto segue, a ritroso rispetto alla tua definizione, ma sono obbligato a sbagliare se provo a cambiare.

  • \phi - l'angolo dall'asse z
  • \theta - l'angolo di proiezione nel piano x-y.

La parametrizzazione

Sia t un insieme discreto di N punti equidistanti da 0 a pi (compreso).

\phi(t) = t 
\theta = 2 * c * t 

Piuttosto dritto in avanti e semplice, una spirale attorno ad una sfera è lineare nel \phi e theta.c è una costante che rappresenta il numero di rotazioni complete in \theta, non deve essere un numero intero.

punti vicini

nel tuo esempio, si calcola l'angolo tra i vettori con atan2(norm(cross....) che va bene, ma non fornisce alcuna comprensione del problema. Il tuo problema è sulla superficie di una sfera, usa questo fatto. Quindi considero la distanza tra i punti using this formula

Ora trovi i punti vicini, questi si verificano a t +- dt e theta +- 2pi per qualsiasi cosa accada.

Nel primo caso t +- dt, è facile da calcolare cos(gamma) = 1 - 2 c^2 sin^2(t) dt^2. La dipendenza sin^2(t) è il motivo per cui i poli sono più densi. Idealmente si desidera scegliere theta(t) e phi(t) tali che dtheta^2 * sin^2(phi) sia costante e minimo per soddisfare questo caso.

Il secondo caso è un po 'più difficile e fa apparire i miei commenti su "sconcertare" i tuoi punti. Se scegliamo una N tale che dtheta non divida equamente 2pi, quindi dopo una rotazione completa intorno alla sfera in theta non riesco a finire direttamente sotto un punto precedente. Per confrontare la distanza tra i punti in questo caso, utilizzare delta t in modo che c delta t = 1. Quindi hai delta phi = delta t e delta theta = 2 c delta t - 2pi. A seconda della scelta di c, delta phi potrebbe essere o meno piccolo abbastanza da utilizzare le approssimazioni di piccolo angolo.

Note finali

Dovrebbe essere evidente che c=0 è una linea retta verso il basso la sfera. Aumentando lo c si aumenta la "densità della spirale" ottenendo una maggiore copertura. Tuttavia, si aumenta anche la distanza tra i punti vicini. Dovrai scegliere uno c per il N scelto che rende le due formule di distanza sopra approssimativamente uguali.

EDIT spostato alcune cose da mathbin per la pulizia

+0

Ho pubblicato una versione modificata del mio codice originale - pensi che incorpori i tuoi suggerimenti, o c'è qualcosa che ho trascurato? Per compilare in Matlab/Octave puoi usare http://www.compileonline.com/execute_matlab_online.php –

+0

Ho anche aggiunto la formula che hai suggerito per la distanza tra due punti. –

0
clear all 
close all 

u = pi/2:(-pi/36):-pi/2; 
v = 0:pi/36:2*pi; 

nv = length(v); 
nu = length(u); 

f = myfigure(1); 
ax = myaxes(f,1,1); 

hold on 

for aa = 1:nv 
    tv = v(aa); 
    for bb = 1:nu 
     tu = u(bb); 
     x = sin(tu)*cos(tv); 
     y = cos(tu)*cos(tv); 
     z = sin(tv); 
     plot3(x,y,z,'*') 
    end 
end 

edit: miafigura e myaxes sono funzioni che ho per la creazione di una figura e un assi

0

ho scritto una versione veloce in C che si comporta molto bene dato un numero fisso di punti. Puoi giocarci allo ideone. Se hai un browser abilitato per WebGL (Chrome, Firefox) puoi incollare quei risultati here per vederli tracciati. I poli sono un po 'spenti a causa di alcune approssimazioni integrali usate nel derivare le formule, ma a parte questo è difficile vedere un difetto. Non ci sono costanti che dovrebbero necessitare di un pizzicamento diverso dal numero di punti per cui si desidera produrre.

#include <stdio.h> 
#include <math.h> 

int main(void) { 
    int i, numPoints = 200; 
    double slope = sqrt(1.2)/sqrt(numPoints); 
    for (i = 0; i < numPoints; i++) { 
     double s = asin((double)i/(double)(numPoints - 1) * 2.0 - 1.0); 
     double z = sin(s); 
     double r = sqrt(1.0 - z * z); 
     double ss = (2.0 * s + M_PI)/slope; 
     double x = cos(ss) * r; 
     double y = sin(ss) * r; 
     printf("%lf,%lf,%lf,\"1\"\n", x, y, z); 
    } 
    return 0; 
}