Per determinare se un punto si trova all'interno di un cubo si può verificare se:
[*] la coordinata "
x" del
punto è
compresa tra la coordinata
"x" minima e la coordinata
"x" massima del
cubo
[*] la coordinata "
y" del
punto è
compresa tra la coordinata
"y" minima e la coordinata
"y" massima del
cubo
[*] la coordinata "
z" del
punto è
compresa tra la coordinata
"z" minima e la coordinata
"z" massima del
cubo
Se le
tre condizioni sono
verificate contemporaneamente, il
punto si trova all'
interno (o su un vertice, spigolo, faccia) del cubo
Si può quindi impostare uno script che:
[*] identifichi il lato dei cubi: deve essere un sottomultiplo degli intervalli nei quali sono comprese le coordinate dei punti (nel caso della domanda (X [500 3000], Y [-100 250], Z [800 1000]).
[*] scorra tutti i punti in loop "esterno"
[*] scorra tutti i cubi in un loop "interno" (al precedente)
[*] verifichi la condizione sopra indicata
[*] registri quali cubi risultano contenere dei punti
Lo "scorrimento" dei cubi, nel loop interno, può essere effettuato incrementando i valori delle coordinate minima e massima in X, Y e Z con un passo uguale al lato scelto inizialmente per il cubo.
La soluzione proposta è stata implementata nello script che segue.
%
% Apertura della prima finestra che conterrà solo i punti
%
figure;
%
% Definizione limiti
%
% X=[500 3000];
% Y=[-100 250];
% Z=[800 1000];
X=[0 400];
Y=[0 400];
Z=[0 400];
dim_x=(X(2)-X(1));
dim_y=(Y(2)-Y(1));
dim_z=(Z(2)-Z(1));
sort_dim=sort([dim_x dim_y dim_z]);
passo=min([gcd(sort_dim(1),sort_dim(2)) gcd(sort_dim(1),sort_dim(3)) gcd(sort_dim(2),sort_dim(3))]);
%
% Aggiustamento del passo
%
passo_manuale=100;
if(passo_manuale)
passo=passo_manuale;
else
fattore_rid=2;
if(find(sort_dim == passo))
passo=passo/fattore_rid;
end
end
%
% Definzione numero punti
%
n_pt=33;
x=randi(dim_x,1,n_pt)+X(1);
y=randi(dim_y,1,n_pt)+Y(1);
z=randi(dim_z,1,n_pt)+Z(1);
%
% Creazione assi cartesiani
%
axes
set(gca,'xlim',X,'ylim',Y,'zlim',Z)
% l_off=.3;
% set(gca,'xlim',[X(1)-floor((abs(X(1))*l_off)) X(2)+floor((abs(X(2))*l_off))], ...
% 'ylim',[Y(1)-floor((abs(Y(1))*l_off)) Y(2)+floor((abs(Y(2))*l_off))], ...
% 'zlim',[min([0,Z(1)-floor((abs(Z(1))*l_off))]) Z(2)+floor((abs(Z(2))*l_off))]);
grid on
view([33 33])
xlabel('Asse X')
ylabel('Asse Y')
zlabel('Asse Z')
hold on
%
% Plot dei punti
%
for i=1:n_pt
p(i)=plot3(x(i),y(i),z(i),'o','markerfacecolor','r');
end
for i=1:n_pt
set(p(i),'visible','on');
end
%
% Calcolo lati, numero lati e numero sezioni per l'asse x
%
lato_x=X(1):passo:X(2);
n_lati_x=length(lato_x);
sez_x=n_lati_x-1;
%
% Calcolo lati, numero lati e numero sezioni per l'asse y
%
lato_y=Y(1):passo:Y(2);
n_lati_y=length(lato_y);
sez_y=n_lati_y-1;
%
% Calcolo lati, numero lati e numero sezioni per l'asse z
%
lato_z=Z(1):passo:Z(2);
n_lati_z=length(lato_z);
sez_z=n_lati_z-1;
%
% Calcolo numero totale (massimo) di cubi
%
n_cubi=sez_x*sez_y*sez_z;
%
% Inizializzazioine array di strutture "cubo"
% Campi:
% pos=vettore indici posizione
% n_pt=numero punti
% pc_1_pt=percentuale (0-1) di punti all'interno del sigolo cubo
% pc_100_pt=percentuale (0-100) di punti all'interno del sigolo cubo
% pc_100_r_pt=percentuale arrotondata (0-100) di punti all'interno del sigolo cubo
% pox_x=indice "x" del cubo
% pox_y=indice "y" del cubo
% pox_z=indice "z" del cubo
% p_x=vettore coordinate "x" dei punti all'interno del singolo cubo
% p_y=vettore coordinate "y" dei punti all'interno del singolo cubo
% p_z=vettore coordinate "z" dei punti all'interno del singolo cubo
%
cubo(n_cubi,1).n_pt=0;
for i=1:n_cubi
cubo(i).pos=[NaN NaN NaN];
cubo(i).n_pt=0;
cubo(i).pc_1_pt=0;
cubo(i).pc_100_pt=0;
cubo(i).pc_100_r_pt=0;
cubo(i).pos_x=0;
cubo(i).pos_y=0;
cubo(i).pos_z=0;
cubo(i).p_x=0;
cubo(i).p_y=0;
cubo(i).p_z=0;
end
%
% Inizializzazione calocolo tempo di esecuzione
%
tic
%
% Identificazione del cubo all'interno del quale si trova ogni punto
% Loop esterno: punti
% Loop interni: cubi
%
for i=1:n_pt
set(p(i),'markerfacecolor','g','visible','on')
for j=1:n_lati_x-1
if(x(i) >= lato_x(j) && x(i) <= lato_x(j+1))
pos_x=j;
break;
end
end
for j=1:n_lati_y-1
if(y(i) >= lato_y(j) && y(i) <= lato_y(j+1))
pos_y=j;
break;
end
end
for j=1:n_lati_z-1
if(z(i) >= lato_z(j) && z(i) <= lato_z(j+1))
pos_z=j;
break;
end
end
%
% Assegnazione output ai campi della struttura "cubo"
%
cubo_idx=(pos_y-1)*sez_x+pos_x+(pos_z-1)*sez_x*sez_y;
cubo(cubo_idx).pos=[pos_x pos_y pos_z];
cubo(cubo_idx).n_pt=cubo(cubo_idx).n_pt+1;
cubo(cubo_idx).pos_x=pos_x;
cubo(cubo_idx).pos_y=pos_y;
cubo(cubo_idx).pos_z=pos_z;
cubo(cubo_idx).p_x(cubo(cubo_idx).n_pt)=x(i);
cubo(cubo_idx).p_y(cubo(cubo_idx).n_pt)=y(i);
cubo(cubo_idx).p_z(cubo(cubo_idx).n_pt)=z(i);
%
% Cambiamento colore del marker del punto (solo per debug)
%
set(p(i),'markerfacecolor','k','visible','on')
end
%
% Calcolo del tempo di esecuzione
%
toc
%
% Calcolo percentuale di punti in ogni cubo
%
for i=1:n_cubi
cubo(i).pc_1_pt=cubo(i).n_pt/n_pt;
cubo(i).pc_100_pt=cubo(i).n_pt/n_pt*100;
cubo(i).pc_100_r_pt=round(cubo(i).n_pt/n_pt*100);
end
%
% Limiti degli assi della prima figura, da assegnare alla seconda figura
%
x_lim=get(gca,'xlim');
y_lim=get(gca,'ylim');
z_lim=get(gca,'zlim');
%
% Apertura della seconda figura che conterrà i punti ed i cubi all'interno dei quali si trovano
%
figure
%
% Creazione assi cartesiani
%
axes
xlabel('Asse X')
ylabel('Asse Y')
zlabel('Asse Z')
grid on
view([33 33])
hold on
set(gca,'xlim',x_lim,'ylim',y_lim,'zlim',z_lim)
%
% Definizione
% trasparenza delle facce dei cubi
% colore delle facce dei cubi
% colore dei markers dei punti
%
trasparenza=.15;
colore='g';
pt_col='k';
%
% Plot dei cubi all'interno dei quali si trova almeno un punto
%
for i=1:n_cubi
if(cubo(i).n_pt)
%
% Display delle informazioni relatice al cubo "corrente" (solo per debug)
%
cubo(i);
x1=X(1)+(cubo(i).pos_x-1)*passo;
x2=X(1)+(cubo(i).pos_x)*passo;
y1=Y(1)+(cubo(i).pos_y-1)*passo;
y2=Y(1)+(cubo(i).pos_y)*passo;
z1=Z(1)+(cubo(i).pos_z-1)*passo;
z2=Z(1)+(cubo(i).pos_z)*passo;
patch('xdata',[x1 x2 x2 x1],'ydata',[y1 y1 y2 y2], ...
'zdata',[z1 z1 z1 z1],'facecolor',colore,'facealpha',trasparenza)
patch('xdata',[x1 x2 x2 x1],'ydata',[y1 y1 y2 y2], ...
'zdata',[z2 z2 z2 z2],'facecolor',colore,'facealpha',trasparenza)
patch('xdata',[x1 x2 x2 x1],'ydata',[y1 y1 y1 y1], ...
'zdata',[z1 z1 z2 z2],'facecolor',colore,'facealpha',trasparenza)
patch('xdata',[x1 x2 x2 x1],'ydata',[y2 y2 y2 y2], ...
'zdata',[z1 z1 z2 z2],'facecolor',colore,'facealpha',trasparenza)
patch('xdata',[x1 x1 x1 x1],'ydata',[y1 y2 y2 y1], ...
'zdata',[z1 z1 z2 z2],'facecolor',colore,'facealpha',trasparenza)
patch('xdata',[x2 x2 x2 x2],'ydata',[y1 y2 y2 y1], ...
'zdata',[z1 z1 z2 z2],'facecolor',colore,'facealpha',trasparenza)
%
% Plot dei punti all'interno del cubo "corrente"
%
n_punti=cubo(i).n_pt;
for j=1:n_punti
p(i)=plot3(cubo(i).p_x(j),cubo(i).p_y(j),cubo(i).p_z(j), ...
'o','markeredgecolor',pt_col,'markerfacecolor',pt_col);
end
end
end
%
% Grafici riassuntivi
%
x=1:sez_x;
y=1:sez_y;
figure
[X,Y]=meshgrid(x,y);
xy=sez_x*sez_y;
for i=1:sez_z
% a1=[cubo(1+xy*(i-1):xy*i).n_pt];
a1=[cubo(1+xy*(i-1):xy*i).pc_100_pt];
a11=reshape(a1,sez_x,sez_y)';
surf(X,Y,ones(size(X))+lato_z(i),a11)
shading interp
hold on
end
caxis([0 max([cubo(:).pc_100_pt])])
% caxis([0 100])
xlabel('Asse X')
ylabel('Asse Y')
zlabel('Asse Z')
title('Percentuale (0-100) Punti per Cubo (raggruppati per alteza)')
colorbar
Hope this helps.
Allegati:
Punti di eempio
Rappresentazione grafica dei cubi contenenti almeno un punto
Distribuzione percentuale dei punti raggruppati per altezza