MATLAB Tricks：提取论文曲线图中的数据（II）

在上一篇博客中，我们首先讨论了论文以及产品手册中常见的曲线图的大致分类，并分析了每类图片的特点。以此为基准，我们将这些图片分为两类处理。对于第一类图片，由于曲线变化较为简单，在曲线上所需提取的数据点数量较少，我们讨论了一种简单直接的数据提取办法，并展示了其程序示例。
而对于变化趋势更复杂的曲线，或者数据点数量要求高的曲线；我们则需要考虑一种更为自动化的数据提取方式。
本篇博客，我们将首先从上文所述的第iii类图片——复杂单曲线图出发，讨论计算机编程自动复原曲线的基本方法；再从此拓展至第iv类图片——复杂多色曲线图，给出适用于上述两类图片的曲线图数据提取完整程序。

2.1 复杂单曲线图的数据提取

i) 程序的设计思路

对于复杂的曲线图，手动选点的方法会花费相当多的时间；并且随着重复劳动的增加，精度也难以保证。因此，我们希望，能够用一种更为自动化的方式，而非上述第一种方法，来获得复原的数据。除此之外，对于第ii类曲线图，我们采用这样的新方法进行数据的提取，也能极大的提高我们工作和研究的效率。

图1-1 单曲线图示例

我们即以如上的图1-1作为我们处理该类问题的示例考虑。首先定义图像复原所需的一系列参数，并使用imread()函数读取图片。MATLAB在读取常见的 .jpg文件 会将图片转化成一个n×m×3的矩阵处理，其中第一、第二维分别对应每一个像素的位置编码，第三维则对应RGB三通道的颜色灰度值。对于单曲线问题，我们无需使用三通道颜色既可以对目标曲线和图中的其他要素进行区分。因此，可以使用rgb2gray()函数将曲线图转变为灰度图进行处理，提升程序的效率，降低编程的复杂度。
此时，获得的变量img_g为一n×m的二维矩阵，其每个单元的数值为unit8格式，变化范围为0至255。在对两个unit8格式的数a与b求差时，若a>b，则a-b>0；若a≤b，则a-b=0。这一格式会使后续的求差问题变得复杂，因此需要进一步转换格式来方便后续处理，此处使用了im2double()函数，将图片数值格式由unit8格式转换为double格式，转换后，像素值的变化范围为0至1。

%程序初始化
clear all;
%定义图像处理所需的参数容限
Tol_grayscale=15/255;                           %输入灰度容差
Tol_LineW=5;                                    %输入线宽容差
%定义真实坐标轴
Lx_real=7;                                      %输入X轴真实坐标长度；
Ly_real=3.5*10^7;                               %输入Y轴真实坐标长度；
x0_real=1;                                      %输入X轴真实起始坐标；
y0_real=0;                                      %输入Y轴真实起始坐标；img_1=imread('Pic_Test1.jpg');                   %读取需要提取数据的图片
img_g=rgb2gray(img_1);
img_g=im2double(img_g);
imshow(img_g,'InitialMagnification',200);       %显示读取的图片，并设置显示比例

随后，我们使用ginput()函数在图片上选取曲线图的四顶点，以标定图片中坐标轴的位置。
其次，我们需要使用ginput()函数在目标曲线上采样，提取一个基准灰度值，以此为基准完成目标曲线上所有像素点位置的提取。对于灰度图片的像素点，灰度值为255时，像素点为纯白色，灰度值为0时，像素点为纯黑色。在double格式下，数值为1时，像素点为纯白色，数值为0时，像素点为纯黑色。一般来说，对于一条曲线上的像素点：靠近边缘部分的像素点颜色较浅、灰度值更高，而曲线中心位置的像素点颜色较深、灰度值更浅。因此，我们尝试在曲线上提取五个点，并对这五个点进行一定的运算，以获得一个更为准确的曲线灰度基准值。使用ginput()函数提取灰度基准采样点的过程中，可能出现误采样而提取非曲线上像素点的可能性，因此不能直接对五个点取平均值以获得灰度基准。我们利用mean()函数计算这五个点的灰度平均值，并利用std()函数计算五点的标准差；剔除mean()±std()以外的值，随后再次计算平均值以获得更为准确的灰度基准值。

%提取坐标轴
disp('选取坐标图的4顶点...');
[x_ax,y_ax] = ginput(4);                        %使用ginput函数实现坐标轴4端点（左下、左上、右下、右上）的图上坐标提取，在图窗上依次点击4端点以完成提取
x_ax=sort(x_ax);                                %使用sort函数对提取的4个x坐标进行排序，从而区分x轴的起点与x轴的终点，前两个数据为x轴的起点，后两个数据为x轴的终点
y_ax=sort(y_ax);                                %同上
x0_pic=round(mean(x_ax(1:2)));                  %计算图像中的X轴起点；
x1_pic=round(mean(x_ax(3:4)));                  %计算图像中的X轴终点；
Lx_pic=x1_pic-x0_pic;                           %计算图像中的X轴长度；
y0_pic=round(mean(y_ax(1:2)));                  %计算图像中的Y轴起点；
y1_pic=round(mean(y_ax(3:4)));                  %计算图像中的Y轴终点；
Ly_pic=y1_pic-y0_pic;                           %计算图像中的Y轴长度；%提取数据
%提取曲线灰度
disp('在目标曲线上选取5个点，获得目标曲线的基准灰度值...');
[x_cr,y_cr]=ginput(5);                  %提取曲线上5个点，作为基准灰度参考值；
gray_cr=zeros(1,5);                     %获取上述五点的灰度值；
for i=1:5gray_cr(i)=img_g(round(y_cr(i)),round(x_cr(i)) );
end
gray_mean=mean(gray_cr);                %计算上述五点的灰度平均值；
gray_std=std(gray_cr);                  %计算上述五点的灰度标准差；
gray_rf=0; k=1;                         %初始化参考灰度值；取标准差以内的灰度值计算平均数获得参考灰度值；
for i=1:5if abs(gray_cr(i)-gray_mean)<=gray_stdgray_rf=(gray_rf*(k-1)+gray_cr(i))/k;k=k+1;end
end

这之后，我们以上面得到的这一灰度基准值为准，以手动设置的参数Tol_grayscale为误差容限，使用find()函数在图像中坐标轴圈定的范围内寻找所有灰度值与该值相近的点。find()函数会返回上述点的x与y坐标。由于图片中的曲线具有一定的宽度，并且坐标轴圈定的图像部份内可能出现一些灰度相近的噪点（例如坐标轴的边缘或者刻度线）；一个x值可能对应多个y点。故而，我们需要用一定的方法，对这多个y值进行辨识，以获得每一个x坐标所对应的最准确的y值，如下：

%在图片中寻找灰度相近的点并记录
[yc_pic,xc_pic]=find(abs(img_g(y0_pic:y1_pic,x0_pic:x1_pic)-gray_rf)<=Tol_grayscale);

一般的，对这一组多个y值，如果噪点数量较小，经常用计算基准灰度类似的方法、取平均数和标准差以剔除两侧因多余框线而出现的噪点。但在实际操作过程中，两侧框线所提供的y值数量甚至与目标曲线提供的y值数量相当，这样，计算得到的初步标准差将会非常大，无法进行有效的误差剔除。因此，我们需要更有效的算法，来获得更准确的y值。我们选择每一组数据的中位数来作为y值的初步基准；再根据手动设定的线宽容限Tol_LineW（这一值一般取实际线宽的两到三倍）来剔除两侧框线带来的多余数据；最后，将筛选得到的数据取平均，以获得最终的y值，如下：

%对find()函数获得的点进行筛选合并；
x_curve=zeros(1,Lx_pic);
y_curve=zeros(1,Lx_pic);
count1=1;
for i=2:length(xc_pic)if (i==2) || (xc_pic(i)~=xc_pic(i-1))           %将x(i)相等的点的y(i)值置于同一数组内x_temp=xc_pic(i);y_temp=yc_pic(i);for j=i+1:length(xc_pic)if xc_pic(i)==xc_pic(j)y_temp=[y_temp,yc_pic(j)];elsebreak;end        endyt_median=median(y_temp);                   %计算y(i)数组的中位数yt_0=0;n2=1;for k=1:length(y_temp)                      %剔除容限外的y(i)值if abs(y_temp(k)-yt_median)<=Tol_LineWyt_0=(yt_0*(n2-1)+y_temp(k))/n2;n2=n2+1;end        end        x_curve(count1)=x_temp;y_curve(count1)=yt_0;count1=count1+1;  end
end

图1-2 基于中位数与线宽容限复原得到的曲线图

即使采用上述方法对噪点进行剔除，实际获得的曲线依然会有个别突变的噪点，如图8所示。这些噪点一般因坐标轴上刻度线凸起的影响而产生。因此，我们还需要在上述过滤算法得到的曲线上进一步滤除突变点，以获得更为准确的复原数据与曲线。突变点的过滤仅需设置一个突变门限，并检测曲线前后两点间y坐标的差距是否超过这一门限即可实现。

%剔除曲线中的突变点
xc_smooth=0;
yc_smooth=Ly_pic;
for i=1:length(x_curve)if abs(y_curve(i)-yc_smooth(end))<50xc_smooth=[xc_smooth,x_curve(i)];yc_smooth=[yc_smooth,y_curve(i)];end
end

最后，我们使用一组线性坐标变换，既可以完成数据的复原。

%剔除曲线中的突变点
xc_smooth=0;
yc_smooth=Ly_pic;
for i=1:length(x_curve)if abs(y_curve(i)-yc_smooth(end))<50xc_smooth=[xc_smooth,x_curve(i)];yc_smooth=[yc_smooth,y_curve(i)];end
end

然而，这样的数据复原仍然存在问题：在数据曲线的起始部分，由于左侧y坐标轴上包含多个刻度凸起，这一部分的曲线坐标识别相当困难，往往最终计算得到的结果仅仅为坐标轴的中点或上下极点坐标，如图1-3所示。

图1-3 左侧刻度突起造成的复原误差

对于这一问题，最简单的处理办法则是，在程序的开始，手动标定这一坐标轴区域的几何位置，并将其灰度值置为1，以消除干扰。代码块如下：

%消除左侧坐标轴突起的干扰
disp('为消除左侧坐标轴突出标识的干扰，选择左侧坐标轴突出标识区的4顶点...');
[x_Lax,y_Lax]=ginput(4);
x_Lax=sort(x_Lax);
y_Lax=sort(y_Lax);
x0_L=round(mean(x_Lax(1:2)));
x1_L=round(mean(x_Lax(3:4)));
y0_L=round(mean(y_Lax(1:2)));
y1_L=round(mean(y_Lax(3:4)));
img_g(y0_L:y1_L,x0_L:x1_L)=1;imshow(img_g,'InitialMagnification',150);

ii) 程序的运行结果

以上，即对程序设计的基本思路和代码作用进行了分步讲解，下图即展示了按上述思路编程得到的数据复原情况：

图1-4 原曲线图与复原结果间的对比

复原得到的数据，只要对数据的头与尾进行简单的处理即可；这样，大大加快了数据复原的效率。

2.2 复杂多色曲线图的数据提取

而对于复杂多色曲线，我们则可以拓展修改上一节所展示的方法：将彩色曲线图按RGB分解成三张灰度图分别进行处理，并使用find()函数在图片中寻找在RGB三个色彩分量上均满足灰度容差要求的点，以此实现对特定颜色曲线的数据提取与复原。

我们以如下的曲线图展示这一数据提取程序的运行效果。在下图中，包含蓝色及橙色两条曲线。

图2-1 彩色多曲线图复原示例

程序对其中蓝色曲线的复原结果如下：

图2-2 蓝色曲线的复原结果

图2-3 橙色曲线的复原结果

可以看到，如上的编程思路可以有效实现包含多条曲线的彩色图片的数据复原，达到了我们预想的要求。

这一程序的完整示例已作为资源上传，链接如下：
https://download.csdn.net/download/weixin_43723517/13138607

以上就是MATLAB Tricks专栏关于论文曲线图复原的全部内容，欢迎读者们交流讨论~