Matlab:“X 射线"通过补丁绘制线 [英] Matlab: "X-Ray" plot line through patch

问题描述

问题

我正在尝试可视化 3-D 路径，以及它周围的云"，它代表数据的标准偏差.我希望能够看到一条粗黑线作为路径，它周围有一个均匀的灰色区域，线条没有任何混浊，就像通过云层看到的线条像 X 射线一样.

I'm attempting to visualize a 3-D path, along with a "cloud" around it which represents the standard deviation of data. I would like to be able to see a thick black line as the path with an EVENLY gray area around it, without any cloudiness to the lines, as if seeing the line through the cloud like an x-ray.

尝试

我使用 plot3 创建一条粗线，使用 patch 创建一系列以线的每个点为中心的框(图中还有一些额外的东西代表开始/停止和方向，我也希望它们很容易被看到).我尝试使用 patch 的 alpha，但这会导致线条混浊，因此灰色的亮度会随着其中有多少灰色框而变化视线.我希望 alpha 为 1，以便每个灰色框的颜色完全相同，但我希望找到某种方法使线条均匀地穿过云层.

I used plot3 to create a thick line and patch to create a series of boxes centered around each point of the line (there are some extra things on the plot to represent start/stop and direction, I would also like them to be seen easily). I tried playing with the alpha of the patch, but that creates a cloudiness to the line, such that the brightness of the gray color changes with however many gray boxes are in the line of sight. I would like the alpha to be 1, so that each gray box is exactly the same color, but I was hoping to find some way to make line seen through the cloud evenly.

最小示例

根据要求，这里是一个最小的例子，它产生下面的图.

As requested, here is a minimal example, which produces the plot below.

% Create a path as an example (a circle in the x-y plane, with sinusoidal deviations in the z-axis)
t = 0:1/100:2*pi;
x = sin(t);y = cos(t);
z = cos(t).*sin(5*t);
figure;
plot3(x,y,z,'k','linewidth',7);

% Draw patches
cloud = .1*rand(size(t)); % The size of each box (make them random, "like" real data)
grayIntensity = .9; % Color of patch
faceAlpha = .15; % Alpha of patch

for i = 1:length(x)
patch([x(i) - cloud(i); x(i) + cloud(i); x(i) - cloud(i); x(i) + cloud(i); x(i) - cloud(i); x(i) + cloud(i); x(i) - cloud(i); x(i) + cloud(i)],... % X values
[y(i) - cloud(i); y(i) - cloud(i); y(i) + cloud(i); y(i) + cloud(i); y(i) - cloud(i); y(i) - cloud(i); y(i) + cloud(i); y(i) + cloud(i)],... % Y values
[z(i) + cloud(i); z(i) + cloud(i); z(i) + cloud(i); z(i) + cloud(i); z(i) - cloud(i); z(i) - cloud(i); z(i) - cloud(i); z(i) - cloud(i)],... % Z values
grayIntensity*ones(1,3),... % Color of patch
'faces', [1 2 4 3;5 6 8 7;1 2 6 5; 8 7 3 4;1 5 7 3;2 6 8 4],... % Connect vertices to form faces (a box)
'edgealpha',0,... % Make edges invisible (to get continuous cloud effect)
'facealpha',faceAlpha); % Set alpha of faces
end

抱歉 for 循环中的代码很长，patch 命令有相当多的参数.前三行简单地定义了定义立方体的 8 个顶点的 x、y 和 z 坐标，通过指定中心点加上或减去立方体的某个半宽，cloud(i).其余的应该由他们各自的评论来解释.

Apologies for the VERY long stretch of code in the for loop, there are quite a large number of arguments to the patch command. The first three lines are simply defining the x, y, and z coordinates of the 8 vertices which define a cube, by specifying the center point plus or minus some half-width of the cube, cloud(i). The rest should be explained by their respective comments.

感谢您的帮助！

推荐答案

这是我在评论中提到的解决方案的实现(只绘制一个 surface 云).

Here is an implementation of the solution I was mentioning in the comment (just drawing one single surface cloud).

它没有优化，有一些 for 循环可以通过巧妙地使用 bsxfun 或这些辅助函数系列来避免，但它运行正常.在每个点处找到曲线的切线和定向(旋转)每个横截面的数学计算也可能会被简化，但这不是我的强项，所以如果他们愿意，我会把它留给专家.

It is not optimized, there are a few for loop which could be avoided with clever use of bsxfun or these family of helper function, but it runs ok as it is. The math for finding the tangent of the curve at each point and orienting (rotating) each cross section could probably be simplified as well, but this is not my strong suit so I leave that to the experts if they feel like it.

基本上，它定义一个圆(在代码中通常称为横截面")，其半径与某物(应用程序中的标准偏差，示例中的随机值)成正比.然后在 3D 中旋转每个圆，使其在平移点处与曲线垂直.然后将所有这些圆用作单个 surface 图形对象的包络.

Basically, it define a circle (often referred as "cross section" in the code) with a radius proportional to something (the standard deviation in your application, a random value in the example). Each circle is then rotated in 3D so it is normal to the curve at the point where it is translated. All these circles are then used as an envelope for a single surface graphic object.

当表面多次重叠时(取决于视角)，主中心线仍有一点阴影，但主线始终可见.此外，您只需管理一个图形对象.

There is still a bit of shadowing of the main centre line when the surface overlap itself many time (depending on the view angle), but the main line is always visible. Plus you only have one graphic object to manage.

结果如下:

当然，您可以根据自己的喜好改变表面的 AlphaValue.我定义了一个与颜色信息数据大小相同的完整矩阵.目前全部设置为0(所以它指向默认colormap中的绿色)，但是这种方式如果你想制作另一个参数的颜色函数也很容易，只需相应地调整颜色矩阵(以及随之而来的颜色图).

Of course you can vary the AlphaValue of the surface to your liking. I defined a full matrix the same size as the data for the colour information. At the moment it is all set to 0 (so it points to the green colour in the default colormap), but this way it is also easy if you want to make the colour function of another parameter, just adjust the colour matrix accordingly (and the colormap that goes with it).

代码末尾有一个选项可以将每个横截面显示为补丁对象.它不打算在最终结果中使用，但如果您想进行自己的修改，它可以帮助您了解整个事物的构造方式.

There is an option at the end of the code to display each cross section as a patch object. It is not intended to be used in the final result, but it can help you understand how the whole thing is constructed if you want to make your own modifications.

代码如下:

%% // Create a path as an example (a circle in the x-y plane, with sinusoidal deviations in the z-axis)
nPts = 180 ;
t = linspace(0,359,nPts)*pi/180;
x = sin(t); y = cos(t);
z = cos(t).*sin(2*t);

figure;
h.line = plot3(x,y,z,'k','linewidth',2,'Marker','none');
hold on
xlabel('X')
ylabel('Y')
zlabel('Z')

%% // Define options
%// cloud = .1*rand(size(t)) ; % The size of each box (make them random, "like" real data)
%// I used another randomization process, make that function of your stdev
r.min = 0.1 ; r.max = 0.2 ;
radius = r.min + (r.max-r.min).* rand(size(t)) ;

%// define surface and patch display options (FaceAlpha etc ...), for later
surfoptions  = {'FaceAlpha',0.2 , 'EdgeColor','none' , 'EdgeAlpha',0.1 , 'DiffuseStrength',1 , 'AmbientStrength',1 } ;
patchoptions = {'FaceAlpha',0.2 , 'EdgeColor','k'    , 'EdgeAlpha',0.2 , 'DiffuseStrength',1 , 'AmbientStrength',1 } ;
patchcol     = [1 0 0] ;  % Color of patch

%% // get the gradient at each point of the curve
Gx = diff([x,x(1)]).' ;                       %'//damn StackOverflow prettifier 
Gy = diff([y,y(1)]).' ;                       %'//damn StackOverflow prettifier 
Gz = diff([z,z(1)]).' ;                       %'//damn StackOverflow prettifier 
%// get the middle gradient between 2 segments (optional, just for better rendering if low number of points)
G = [ (Gx+circshift(Gx,1))./2 (Gy+circshift(Gy,1))./2 (Gz+circshift(Gz,1))./2] ;

%% // get the angles (azimuth, elevation) of each plane normal to the curve

ux = [1 0 0] ;
uy = [0 1 0] ;
uz = [0 0 1] ;

for k = nPts:-1:1 %// running the loop in reverse does automatic preallocation
    a = G(k,:) ./ norm(G(k,:)) ;
    angx(k) =  atan2( norm(cross(a,ux)) , dot(a,ux))  ;
    angy(k) =  atan2( norm(cross(a,uy)) , dot(a,uy))  ;
    angz(k) =  atan2( norm(cross(a,uz)) , dot(a,uz))  ;

    [az(k),el(k)] = cart2sph( a(1) , a(2) , a(3) ) ;
end
%// adjustment to be normal to cross section plane the way the rotation are defined later
az = az + pi/2 ; 
el = pi/2 - el ;

%% // define basic disc
discResolution = 20 ;
tt = linspace( 0 , 2*pi , discResolution ) ;
xd = cos(tt) ;
yd = sin(tt) ;
zd = zeros( size(xd) ) ;

%% // Generate coordinates for each cross section

ccylX = zeros( nPts , discResolution ) ;
ccylY = zeros( nPts , discResolution ) ;
ccylZ = zeros( nPts , discResolution ) ;
ccylC = zeros( nPts , discResolution ) ;

for ip = 1:nPts
    %// cross section coordinates, with radius function of [rand] in this
    %// example. Make it function of the stdev in your application.
    csTemp = [ ( radius(ip) .* xd )  ; ... %// X coordinates
               ( radius(ip) .* yd )  ; ... %// Y coordinates
                               zd    ] ;   %// Z coordinates

    %// rotate the cross section (around X axis, around origin)
    elev = el(ip) ;
    Rmat = [ 1     0           0     ; ...
             0 cos(elev)  -sin(elev) ; ...
             0 sin(elev)   cos(elev) ] ;
    csTemp = Rmat * csTemp ;

    %// do the same again to orient the azimuth (around Z axis)
    azi = az(ip) ;
    Rmat = [ cos(azi)  -sin(azi) 0 ; ...
             sin(azi)   cos(azi) 0 ; ...
               0            0    1 ] ;
    csTemp = Rmat * csTemp ;

    %// translate each cross section where it should be and store in global coordinate vector
    ccylX(ip,:) = csTemp(1,:) + x(ip) ;
    ccylY(ip,:) = csTemp(2,:) + y(ip) ;
    ccylZ(ip,:) = csTemp(3,:) + z(ip) ;
end

%% // Display the full cylinder
hd.cyl = surf( ccylX , ccylY , ccylZ , ccylC ) ;

%// use that if the graphic object already exist but you just want to update your data:
%// set( hd.cyl , 'XData',ccylX , 'YData',ccylY ,'ZData',ccylZ ) 

set( hd.cyl , surfoptions{:} )


%% // this is just to avoid displaying the patches in the next block
%// comment the "return" instruction or just execute next block if you want
%// to see the building cross sections as patches
return 

%% // display patches
hp = zeros(nPts,1) ;
for ip = 1:nPts
   hp(ip) = patch( ccylX(ip,:) , ccylY(ip,:) , ccylZ(ip,:) , patchcol ) ;
   set( hp(ip) , patchoptions{:} )
end

<小时>

这只是一个带有补丁的快速放大视图(代码重新运行，点数较少，否则它很快就会使整个图形变得混乱):

---

And this is just a quick zoomed view with the patches on (code rerun with a lower number of points, otherwise it quickly clutter the whole figure):

这篇关于Matlab:“X 射线"通过补丁绘制线的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！