Matlab中的点数据生成图

---

简单的背景介绍

今天蠢师弟用comsol导出了一个模拟数据，用matlab打开一看是一个数据长度$189739\times 3$的矩阵

好家伙！点进去一看数据格式还可以，第一列是x坐标，第二列是y坐标，第三列是对应的数据值，但是用普通的画图手段无法将图画出。问了一下师弟，原来这货导出的是comsol中网格点上的数据值，且他的模型用的是三角形网格，无法用一个$M\times N$的矩阵进行表示。好嘛，只能用后处理的手段把图画出来。

尝试plot3

matlab里面的plot3函数可以对坐标和数值进行一一对应画图，优点是出图速度快。缺点是画的点/线图，如果细看的话会有间隙；单条plot3指令只能使用一种颜色，改变色彩麻烦。师弟的这种情况只能采取分层着色的方式进行绘图

load('matlab.mat');
cc=colormap(hsv);%着色空间

z_max = max(u_4(:,3));z_min = min(u_4(:,3));
z = u_4(:,3);
for i = 1:64
    u_ind = find(z_min+(z_max-z_min)*(i-1)/64<z & z<z_min+(z_max-z_min)*i/64);%找到对应的着色空间对应的数据的目录
    plot3(u_4(u_ind,1),u_4(u_ind,2),u_4(u_ind,3),'.','Color',cc(i,:)); hold on;%画出对应的颜色
end
view(0,90)

可以使用一张动图来展示这个过程：

最终的结果是这样的：

看着是不是还不错？但是放大以后可以看到：

网格点不是规整的，是三角形网格。这大大影响到了后面的处理过程，且数据点不方便读取与使用，影响接下来的数据分析。

考虑将数据点重新提取整合成一个矩阵

构造矩阵

构造之前先对数据进行分析：

轴	最小值	最大值
X轴（m）	$-3\times 10^{-3}$	$3\times 10^{-3}$
Y轴 （m）	$-2\times 10^{-3}$	$2\times 10^{-3}$
Z轴（数据值）	41.4421	63.4641

希望是把这个数据放到一个$x_{num}\times y_{num}$的矩阵里面，为方便调试先使用$127\times 127$的矩阵进行保存，然后将数据按照坐标放到矩阵中。由于使用的矩阵数据量是小于总数据量的，说明存在多个数据放到同一个矩阵元素中的情况，需要对矩阵元素取平均，于是有了下面这段小程序：

clear all
close all
load('matlab.mat');
x_max = max(u_4(:,1));x_min = min(u_4(:,1));
y_max = max(u_4(:,2));y_min = min(u_4(:,2));
z_max = max(u_4(:,3));z_min = min(u_4(:,3));
x = u_4(:,1);
y = u_4(:,2);
z = u_4(:,3);
%% 初始化
x_num = 128;
y_num = 128;
cx = linspace(x_min,x_max,x_num);
cy = linspace(y_min,y_max,y_num);
c_data = zeros(x_num-1,y_num-1);
for i = 1:x_num-1
    for j = 1:y_num-1
        u_find = find(cx(i)<x & x<cx(i+1) & cy(j)<y & y<cy(j+1));
        c_data(i,j) = mean(z(u_find));
    end
    disp(num2str(i))
end
figure;surf(c_data');view(0,90);shading interp;colormap(hsv)

结果是下图，可以看出来出来的图相较于plot3画出来的细节少了不少，即分辨率降低了，但是这个是可以通过改变存储矩阵大小来改善的

可是这里面有个致命的问题！即用时太长！

步数	时间（s）
125	4.8471
126	4.8808
127	4.9153

没办法，这个二重循环就是逊啦，如果变成$255\times 255$的矩阵，用时可达19.8258s

遍历填充

无优化版遍历填充

考虑到二重循环实际是对二维矩阵元素寻址以及find函数效率低下的问题，需要解决的问题有：1.提升寻址效率，2.降低或者不使用find函数。

矩阵地址存储是按照顺序排布的，但是在构造矩阵中，使用for循环会破坏地址顺序，导致寻址时间变长，而双重循环会进一步加长这个时间，于是可通过降低循环层数并保留地址顺序来降低总时间消耗。其次，find函数实际是一个比较函数，它需要对数据整体进行遍历，多一个find函数相当于多一重循环，因此，不用或者少用find函数（如使用矩阵运算来代替find函数）能提升效率。

根据这一思路，可以使用遍历数据，填充矩阵的方式来提高运行效率。二话不说，先上代码

clear all
% close all
load('matlab.mat');
x_max = max(u_4(:,1));x_min = min(u_4(:,1));
y_max = max(u_4(:,2));y_min = min(u_4(:,2));
z_max = max(u_4(:,3));z_min = min(u_4(:,3));
x = u_4(:,1);
y = u_4(:,2);
z = u_4(:,3);
%%
x_num = 256;
y_num = 256;
x_delt = (x_max-x_min)/(x_num-1);
y_delt = (y_max-y_min)/(y_num-1);
cx = linspace(x_min,x_max,x_num-1);
cy = linspace(y_min,y_max,y_num-1);
c_data = zeros(x_num,y_num);
c_data_num = c_data;
c2 = nan * zeros(x_num,y_num);
x2 = floor((x-x_min)/x_delt)+1;
y2 = floor((y-y_min)/y_delt)+1;
tic
for i = 1:length(z)
   c_data(x2(i),y2(i)) = c_data(x2(i),y2(i))+z(i);
   c_data_num(x2(i),y2(i)) = c_data_num(x2(i),y2(i))+1;
%    if ~mod(i,10000) disp(num2str(i));end
end
toc
gx = c_data./c_data_num;
figure;surf(gx');view(0,90);shading interp;colormap(hsv)

在上面的程序中，使用c_data 和 c_data_num 两个矩阵分别记录累加数据和累加次数，gx是最后生成的矩阵。保存的矩阵大小为$255\times 255$，而运行的总时间却大大减少了。

时间已过 0.009621 秒。

结果图如下，可以看到边界处有一些不平整点，那是在导出的网格边界处无数据，看不得不平整的话后处理一下就行了。

优化后的遍历填充

其实到这里已经优化了2000多倍了，但是事实上还有优化的空间。在上面的无优化版遍历填充中还是使用了一重循环对数据进行了遍历，这里产生了无必要的寻址时间消耗，可以通过matlab矩阵特性进行优化（其实就一句，但是效率极大的提升了）

d_data = nan * c_data;
d_data((x2-1)*x_num+y2) = z; 

>>时间已过 0.003155 秒。

结果图如下：

与优化前相比，区别在于多个数据放到同一个矩阵元素中时，使用一个数据代替对矩阵元素取平均的过程，因此在细节处会略有不同，也算是一个遗憾吧

遗憾解决了，但是效率降低了

tic
d_data((x2-1)*x_num+y2) = d_data((x2-1)*x_num+y2)+z; 
d_data_num((x2-1)*x_num+y2) = d_data_num((x2-1)*x_num+y2)+1; 
toc
>>时间已过 0.008568 秒。