Info

此问题已关闭。 请重新打开它进行编辑或回答。

I would like to draw a 4d isosurfce plot, where the 4th dimension, is the density of each grid voxel, and I would like it to be thee isovalue

1 次查看（过去 30 天）
Guy Nir
Guy Nir 2019-6-18
关闭： MATLAB Answer Bot 2021-8-20
M is a 3d grid of points.
Ro(t) is the density at grid point t.
Ro3D is thee density at each pixel as a 3D matrix that corrwesponds to each voxel.
How would you plot it?
Thank you!
M=DataSRX(find(idx==1),4:6);
StepSize=100; % [nm]
Xgrid=min(M(:,1)):StepSize:max(M(:,1));
Ygrid=min(M(:,2)):StepSize:max(M(:,2));
Zgrid=min(M(:,3)):StepSize:max(M(:,3));
Sigma=StepSize;
XGsize=numel(Xgrid); YGsize=numel(Ygrid); ZGsize=numel(Zgrid);
Z=zeros(1,XGsize*YGsize*ZGsize);
Ind=zeros(size(M,1),1);
Ro=Z;
% Z(1)=Xgrid(1),Ygrid(1),Zgrid(1); Z(2)=Xgrid(2),Ygrid(2),Zgrid(2);
% Z(Xn)=Xgrid(Xn),Ygrid(1),Zgrid(1); Z(Xn+1)=Xgrid(1),Ygrid(2),Zgrid(1);
Norm=1./((Sigma*sqrt(2*pi)^3));
for i=1:size(M,1)
tx=ceil((M(i,1)-min(M(:,1)))/StepSize); ty=ceil((M(i,2)-min(M(:,2)))/StepSize); tz=ceil((M(i,3)-min(M(:,3)))/StepSize);
if ty==0
ty=1;
end
if tx==0
tx=1;
end
if tz==0
tz=1;
end % if tz==0
t=(XGsize+YGsize)*(tz-1)+XGsize*(ty-1)+tx;
Z(t)=Z(t)+1;
Ind(i)=t;
Ro(t)=Norm*exp(((M(i,1)-Xgrid(tx))^2+(M(i,2)-Ygrid(ty))^2+(M(i,3)-Zgrid(tz))^2)/(2*Sigma^2))+Ro(t);
end % for i=1:size(M,1)
RoDense=Ro.*Z;
Ro3D=reshape(RoDense,[XGsize,YGsize,ZGsize]);
%{
Ro3D=zeros(XGsize,YGsize,ZGsize);
for i=1:numel(Ro)
[I,J,K]=ind2sub([XGsize,YGsize,ZGsize],i);
Ro3D(I,J,K)=Ro(i);
end % for i=1:numel(Ro)
%}
Ro3D_round=round(Ro3D,3);

此问题已关闭。

回答（1 个）

Guy Nir
Guy Nir 2019-6-18
To make it more clear I attached the data, and added code to parse the file.
% Plot density
%% Parse file
FileName='/Users/guynir/Dropbox (HMS)/Collaborations/MAMR Lab/Irene/CSV from paper/CSV_6columns/2017-08-03-19-34-25_Location-05_6column.csv';
[Header, DataSRX] = ReadSRXParticlesCSV(FileName);
figure(); scatter3(DataSRX(:,4),DataSRX(:,5),DataSRX(:,6),20,'filled')
axis equal
%% Separate homologs
Ploidy=2;
idx = kmeans(DataSRX(:,4:6),Ploidy);
figure(); scatter3(DataSRX(find(idx==1),4),DataSRX(find(idx==1),5),DataSRX(find(idx==1),6),20,'filled','r')
hold on
scatter3(DataSRX(find(idx==2),4),DataSRX(find(idx==2),5),DataSRX(find(idx==2),6),20,'filled','b')
axis equal
%% Surface plots
M=DataSRX(find(idx==1),4:6);
shp = alphaShape(M);
figure(); plot(shp)
%% Density map
StepSize=100; % [nm]
Xgrid=min(M(:,1)):StepSize:max(M(:,1));
Ygrid=min(M(:,2)):StepSize:max(M(:,2));
Zgrid=min(M(:,3)):StepSize:max(M(:,3));
Sigma=StepSize;
XGsize=numel(Xgrid); YGsize=numel(Ygrid); ZGsize=numel(Zgrid);
Z=zeros(1,XGsize*YGsize*ZGsize);
Ind=zeros(size(M,1),1);
Ro=Z;
% Z(1)=Xgrid(1),Ygrid(1),Zgrid(1); Z(2)=Xgrid(2),Ygrid(2),Zgrid(2);
% Z(Xn)=Xgrid(Xn),Ygrid(1),Zgrid(1); Z(Xn+1)=Xgrid(1),Ygrid(2),Zgrid(1);
Norm=1./((Sigma*sqrt(2*pi)^3));
for i=1:size(M,1)
tx=ceil((M(i,1)-min(M(:,1)))/StepSize); ty=ceil((M(i,2)-min(M(:,2)))/StepSize); tz=ceil((M(i,3)-min(M(:,3)))/StepSize);
if ty==0
ty=1;
end
if tx==0
tx=1;
end
if tz==0
tz=1;
end % if tz==0
t=(XGsize+YGsize)*(tz-1)+XGsize*(ty-1)+tx;
Z(t)=Z(t)+1;
Ind(i)=t;
Ro(t)=Norm*exp(((M(i,1)-Xgrid(tx))^2+(M(i,2)-Ygrid(ty))^2+(M(i,3)-Zgrid(tz))^2)/(2*Sigma^2))+Ro(t);
end % for i=1:size(M,1)
RoDense=Ro.*Z;
Ro3D=reshape(RoDense,[XGsize,YGsize,ZGsize]);
Ro3D_round=round(Ro3D,3);

此问题已关闭。

标签

产品


版本

R2018b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by