EDA365电子论坛网

标题: 利用matlab从图片中提取曲线坐标数据 [打印本页]

作者: piday123    时间: 2021-1-20 18:57
标题: 利用matlab从图片中提取曲线坐标数据
9 M7 \* ]% V1 ?3 G0 L- L
目录0 M( g& ^3 z( z! J  }
0.引言9 Y; M5 w* H. C% y6 |& `6 Y( W
1.思路详解与分析2 }2 U3 @+ I6 @' \" T
2.MATLAB程序( K# q8 t3 I) b1 e% p1 d' Z
% |  g! @. z4 _$ s
0.引言: A( q( a: v6 @/ M4 p- ?; {+ L. M
  在读文献的时,经常遇到这样的情况:文章里提出的方法好有趣啊,好想拿文中用的数据来试试看看能不能得到相近的结果,可是文中只有根据原始数据绘制的曲线图,没有数据。如下图所示。7 b- v1 `0 l$ X/ E) b* G; ~

  R+ h) \: L( o. w/ w  I( V 2 A( p% I8 b/ ^

' K9 V! H2 i& O' ^6 @0 |* z& q; E  此时,如果能从文中把这幅图截取下来,输入到一个函数中去,最后能返回从图片中提取到的曲线的坐标数据,岂不美哉。这便是本文的工作。
3 N  S' |" M' F: T5 G3 G. ^7 w5 o
7 h' n! n; @; A( p3 W1.思路详解与分析
7 k1 m$ N- J- u( B* S+ a6 O4 ]1.1准备待提取数据的曲线图片
0 @- J5 M7 M" G4 {4 T0 T% i; a. l1 ?$ t5 `
将待提取数据的曲线的图片(如.jpg格式图片),利用 imread 输入到 matlab 中。
: h* O/ Q& B0 [; X  Z/ u/ S
+ n0 _! p* s; H! T$ U- O4 I( |1.2曲线图片预处理与数据转换6 J) p& \/ _2 R/ f$ Q) X# h  u

8 R7 |: r; c) e! m曲线图片预处理步骤的主要工作包含如下:
9 s1 p4 x- q2 Z: \+ g' y
/ G8 W% x; s# R$ l2 e: {& ^(1)图像二值化; |. ^) J- w: m6 O  }

6 ?! R" m* ?8 W* M/ s! J将输入图像进行二值化处理,但分割得到的结果并不全为数据,其中可能还包括坐标轴等干扰点需要去除。
- z' v/ X' {% o& b2 a' i# e2 E. l
(2)获取从图片像素到曲线坐标的定标数据, d6 t6 |# a+ E7 r8 ?

% N5 p: W' |5 v- i2 z" r6 M首先,通过ginput()手动从图片中提取到两个像素点,这两个点分别为曲线坐标框的左上角和右下角。8 {* H2 i( k5 m

+ |1 r" J* I4 t9 y4 s, A 4 S/ Z7 e- f; y) t5 C* O' ^

4 Q5 h) d  X/ W" B& f( {% m此时,便获得了曲线在图片上的像素范围
% P8 R' y1 c* \3 f" k4 [0 R. K! ]) X, j( S: N  ~9 Y! n
[x_index_min, x_index_max] & [y_index_min, y_index_max]" F4 t0 `0 y. t( d/ {) |
' S2 P! c% _4 L% n. |, v4 W
然后,手动输入实际曲线的数据坐标范围 [x_min, x_max] & [ymin, y_max_]. 如下所示。( W: W" U) Z* \7 O8 H( ]+ B
) u, T. g1 f  S
% T  Q6 H& m% Z0 s8 d
% W9 h+ L9 k! G
此时,一方面得到了像素坐标,一方面得到了实际坐标。接下来便利用这对数据,将图片中全部的像素坐标转换到实际坐标。. `% h. T. ~4 i3 U; k$ f  G) {" C
! c4 l) T* C$ h. S
最终,得到了由图片提取到的数据散点图,如下:
9 w% p5 N8 m* L& F5 m, V8 j$ M; ]  q4 J" ]2 t, x1 U

  M. z$ R5 m2 z/ R8 X- E可以看到,此时得到的结果,虽然曲线与所需要的相近,但曲线外的部分,如坐标轴框、坐标轴刻度等也被转换成了数据,还需要进一步的处理。
# ]0 Q$ X0 i  e' b1 g; c6 A0 k- u
" k, N' G: R0 k# f1.3数据的进一步处理并得到最终曲线; k( c9 i* k) b
' ~8 s0 Y, t6 j' U
这一步的主要工作是在数据中除去曲线之外的部分(包括坐标轴框、坐标轴刻度等);以及解决一个x数据对应多个y数据的情况。
1 @' _: \( T3 `+ @1 W4 f4 ?0 K# F) q. Y4 Y2 L* p
显然,坐标轴在整幅数据中,均处于边界位置,因此,很容易想到的一种方法是,设定阈值,将距离边界较近的数据直接删除。这里,设定了两个阈值,一个用来限定x方向上的数据,一个用来限定y方向上的数据。比如设定:rate_x=0.08; rate_y=0.05; 意思是阈值设定为曲线最前端8%和最后段8%的数据,曲线最顶端5%和最底端5%的数据。/ w  H0 J! q5 _: |8 ~
& X1 D" S; T- ~0 X7 [4 ~9 W
进一步的,对于提取到的数据图,大多数情况一个x会对应若干y,因为数据是由图片转化而来,而图片的分辨率有限,一个实际数据点会用多个像素来表示。解决此问题的中心思想是将同一个x对应的若干个y取均值,但不能直接求均值,因为还有很多y是噪声(如坐标轴线、由图片噪声带入的干扰点等)需要先去除,在第一个问题中,通过限定y的范围,已经在一定程度上除去了部分干扰,在此基础上,我们求取一个x对应的这组y值的均值mean与标准差std,当某些y值位于[mean-std , mean+std]之外,则认为这些y值波动太大,将它们删去。
" W3 O" p( N" d2 `8 |9 D/ n
9 N* {5 k6 k/ A8 ?1 w到这里,我们就将数据的处理部分基本完成了,我们将处理后的数据再次绘制成曲线,便可以得到如下
4 @" @9 I; ~% z. n1 K4 @
* J! x" `% H9 O/ Q4 m5 L- i, R # _- J! V% \) C2 R/ D

% ]9 E4 h$ A0 N" E* ?: l对比处理之前的数据,由于限定了范围,因此曲线图片中带来的坐标框等内容转化而来的数据已经被删去。, M, l/ h$ y7 _2 {* U, ?
/ H% u! U9 K, b' @
将需要提取坐标的曲线图片,和我们提取并处理后的数据绘制的曲线,放在一起对比如下:
& _2 _3 q1 ~6 i, f* M1 Z0 k: M
* Y( f( R# z; @& E 9 |5 K1 ~: P) t

; ]! s+ D) w: E. Q# X) a可以看到,与原曲线图片相比,提取到的数据曲线相似度能达到较高要求。但进一步观察会发现,右图曲线较左图而言,高频分量有一定的减少(即右图曲线更平滑),原因在于数据处理时,对同一个x对应的这组y值进行了均值处理,则在图像上近似反应为均值滤波,从而使得提取到的数据绘制成的曲线的高频分量被滤除。" y5 N$ G1 A* q, m3 \7 ]0 V

: Z, H6 m3 a) i5 }( i7 x# i最后,将提取到的最终数据,保存起来如下:
7 x9 _2 q' a; Y1 o' w) W" x- ]' d; P; J3 G0 y" z

$ M' v( w: G. ~8 z" a5 T8 }' w9 S" t, m
1.4进一步的讨论——曲线拟合
7 Q% z. }4 N' w  ~8 m
, Y" Q# y1 L5 n  K- Y# M9 v7 i通过对图片中曲线的数据提取,可以得到数值上的答案,这会带来进一步的思考,即能否得到这些数据的解析表达式。很容易想到,利用最小二乘法来拟合这些数据,这便涉及到了曲线的拟合。(插值与拟合可以这么理解:对于数据点集,若均落在曲线上,则该曲线为插值曲线,否则为拟合曲线)+ C0 @& ?& U) u, b0 E; }4 r5 X
+ I+ R& r* S9 n( a' x7 w3 Y
对于一些简单的曲线图片(如下),可以考虑用泰勒级数来近似,即多项式拟合。
- G0 |3 j; W& i/ G3 o5 r. p* W. n
) q  Y- O4 i, s) Z
# l4 `: _% P5 T6 J8 H
: X; v, l0 d& ~9 k; |" {数据提取并拟合的结果展示如下
3 p+ l+ `4 o: q. s
* r, {- m# ?5 _
( K* [  {3 x" Q! u, X
$ v) F( U# J, E7 O$ t同时还能得到拟合多项式的系数
* }0 o* v% }; j# P" w* I' b* L2 d7 w/ }/ X# W2 q) o
2 I. G5 X/ z# V' N+ ^2 E7 Y
$ Q4 w' O, b# x/ t! Y6 H
从而得到该曲线的多项式(这里采用四阶多项式)表达式为:
- D4 A% C& g2 c: E
- |2 t% V  _' k% Z: ~5 g( J# I   j1 I/ w/ D$ H! q' u) K- j3 F

7 N/ x# `* z" d6 U" h& I4 z理论上,泰勒级数可以分解任何函数,但实际上,多项式拟合的次数太高,会出现龙格库塔现象,即摆尾现象。因此,多项式拟合的阶数不易过高,一般低于5阶。对于本文最开始的那幅曲线而言,仅5阶的泰勒级数就显得力不从心了,因此,对于这种存在波动剧烈的函数,可以考虑用傅里叶级数进行拟合,或者如果能提供先验知识,可以直接用先验表达式进行拟合。
  M/ T, Z  c, @; [/ g2 G7 K2 l- j& b. Z: x
在MATLAB中,提供了cftool工具箱,其提供了拟合与插值的GUI,使用非常方便,直接在命令窗输入cftool即可调用,cftool界面如下所示,其具体使用方法不在此赘述。
- |6 k8 T: Q5 h& y2 {* E$ _- r% |
% a! _, k6 h% o0 |- T 7 _& V6 }* z( @( `& [

; t" V* S! Y9 \( a/ V. R2.MATLAB程序
+ g! J+ u2 S. s; zMATLAB源代码如下所示,和以往的风格一样,提供了详细的注释0 W8 S* Q/ F0 D  b" S

- x. M6 `5 y+ h% //提取图片中的曲线数据/ j, A% o3 [: i$ a$ i3 C4 E7 f
clear,clc,close all
2 J" v# E% c# c1 `%% //图片与曲线间的定标/ c( \4 \* U( Q1 v: j% {
im=imread('tu1.jpg'); %//读入图片(替换成需要提取曲线的图片); m' C8 f# ~" A1 Q6 z6 [! W
im=rgb2gray(im); %//灰度变化2 Q3 @# t. B9 o" R6 C
thresh = graythresh(im); %//二值化阈值2 n5 P5 L) Q! g. j  W
im=im2bw(im,thresh); %//二值化
7 i8 T" Z" f0 I6 ]9 Z' J, r2 kset(0,'defaultfigurecolor','w')
( \' P2 E! w' y$ oimshow(im) %//显示图片
- B4 }4 ^, j. Q- G+ r  K, \& ?[y,x]=find(im==0); %//找出图形中的“黑点”的坐标。该坐标是一维数据。
, M$ \1 ?/ m7 Q/ }9 ^# Z; ]y=max(y)-y; %//将屏幕坐标转换为右手系笛卡尔坐标( [( M4 W0 C9 R, B9 w( e( V# s0 b, B
y=fliplr(y); %//fliplr()——左右翻转数组9 T2 `# b  }9 P2 M% e
plot(x,y,'r.','Markersize', 2);. @* E) U6 M  d
disp('请在Figrure中先后点击实际坐标框的两个顶点(左上点和右下点),即A、B两点. ');" y# ^& n5 W: P$ m
[Xx,Yy]=ginput(2); %//Xx,Yy——指实际坐标框的两个顶点3 ~# l+ D# r* f" c" n# E
min_x=input('最小的x值');  %//输入x轴最小值# K( }! P" A# G# N2 u
max_x=input('最大的x值');  %//输入x轴最大值" L* F1 V; ?$ a+ e7 T$ m7 P3 F
min_y=input('最小的y值');  %//输入y轴最小值6 b; i8 ^6 M! u4 t5 [* I
max_y=input('最大的y值');  %//输入y轴最大值
: g$ C, Q$ a- k) ~! O1 Jx=(x-Xx(1))*(max_x-min_x)/(Xx(2)-Xx(1))+min_x;3 y! r  ^' M8 ?$ ^
y=(y-Yy(1))*(min_y-max_y)/(Yy(2)-Yy(1))+max_y;
& ?/ k6 K3 q$ A0 H7 L6 @" a8 splot(x,y,'r.','Markersize', 2);
/ K( F' [9 K" a7 v7 \! qaxis([min_x,max_x,min_y,max_y])  %//根据输入设置坐标范围' g1 _5 a# d( V4 v- w) c( w
title('由原图片得到的未处理散点图'): V- R( p) c+ u" x  E
%% //将散点转换为可用的曲线
9 M7 E. g' M" u  Z1 n%//需处理的问题与解决思路( }- R; ]5 {: v% d5 Z3 i9 ]
%//(1)散点图中可能一个x对应好几个y <---> 保留mean()-std()到mean()+std()之间的y值 并取平均处理( M$ C3 R6 }! c, i5 U) G8 y9 |0 p$ [
%//(2)曲线的最前端和最后段干扰较大 <---> 去掉曲线整体的前(如5%)和后5%" B0 t( L) @  F% ~
%//(3)曲线的最顶端和最底段干扰较大 <---> 去掉曲线整体的上10%和下10%1 S& o! }  \1 m% o. l
+ B' Q1 ^. q% ^
%//参数预设
- }' M* B; d* Z% b7 p+ qrate_x=0.08;  %//曲线的最前端和最后段删除比例
8 |/ t0 r+ o7 d* Brate_y=0.05;  %//曲线的最顶端和最底段删除比例' L7 V0 P8 M* B$ b4 n! \

( V: f, k( \* [$ B8 T2 O- P) u[x_uni,index_x_uni]=unique(x);  %//找出有多少个不同的x坐标
1 H+ _+ B/ {6 P. V' \; I) \
& ?8 \. S' S5 w0 ^- f' Ix_uni(1:floor(length(x_uni)*rate_x))=[];  %//除去前rate_x(如5%)的x坐标5 l4 x; _6 s& @, |0 ?1 U9 P
x_uni(floor(length(x_uni)*(1-rate_x)):end)=[];  %//除去后rate_x的x坐标  ?1 Y% s0 X3 Z/ j$ _1 H) {, b
index_x_uni(1:floor(length(index_x_uni)*rate_x))=[];  %//除去前rate_x的x坐标+ [& j2 n3 N! w8 r+ {
index_x_uni(floor(length(index_x_uni)*(1-rate_x)):end)=[];  %//除去后rate_x的x坐标  U4 B1 Q6 A. c1 O# R' j

8 C: X* f9 L& O; ^; @9 M[mxu,~]=size(x_uni);
- e5 H0 _5 J6 n/ X1 T[mx,~]=size(x);8 q2 Z2 f0 |" |
for ii=1:mxu
, x8 T6 m; x# n, V    if ii==mxu6 ]- v( k8 u0 ?7 u: H
        ytemp=y(index_x_uni(ii):mx);
. v9 _. F$ D4 z5 Q+ y, x7 l8 y    else
- E# s' c$ L- ?) U( l5 T9 l        ytemp=y(index_x_uni(ii):index_x_uni(ii+1));
/ k7 ~3 |' ]% S* l. i$ I$ b, M' ^  N    end
* _3 ?: [  T' f  S8 k# l    % //删除方差过大的异常点! e" q# N# O$ A& z& K- b: x- I
    threshold1=mean(ytemp)-std(ytemp);
9 D& @& E2 q/ s6 u    threshold2=mean(ytemp)+std(ytemp);2 u2 A- {3 J& s
    ytemp(find(ytemp<threshold1))=[];  %//删除同一个x对应的一段y中的异常点; e! ]9 K3 K  m7 s
    ytemp(find(ytemp>threshold2))=[];
! V. B5 b5 ^9 i4 K8 M* h    % //删除距顶端和底端较近的点8 @. ?! S5 k! a* e
    thresholdy=(max_y-min_y)*rate_y;  %//y坐标向阈值
" i( d, @8 P5 c$ n$ h; ?    ytemp(find(ytemp>max_y-thresholdy))=[];  %//删除y轴向距离顶端与底端距离小于rate_y的坐标) o# F5 |/ @9 r6 V6 M! g0 M
    ytemp(find(ytemp<min_y+thresholdy))=[];+ D6 c0 r' u  j4 g7 ~
    % //剩下的y求均值% \+ R6 x" g2 r; g4 N9 _/ S
    y_uni(ii)=mean(ytemp);# e  s: {2 r' E  y3 o2 V  K
end
$ h6 D' A/ o+ ?# f% //此时很多x_uni点处对应的y_uni为空,即NAN,要进一步删去这些空点4 ^1 I+ Q: F- ?& E7 D4 s
x_uni(find(isnan(y_uni)))=[];7 D! }/ `. j& q3 h) m! M. i' p
y_uni(find(isnan(y_uni)))=[];
3 c7 ?. A: U/ l6 @9 ]& g+ \% //画图# p% e8 j! @+ m+ i1 s0 @7 S
figure,plot(x_uni,y_uni),title('经处理后得到的扫描曲线')8 h9 ]* {$ Z7 k
axis([min_x,max_x,min_y,max_y])  %//根据输入设置坐标范围& @' f  ~5 I1 m4 E/ u2 j
% //将最终提取到的x与y数据保存
% ~# _" s7 `' R( E0 D1 zcurve_val(1,:)=x_uni';5 m+ L, m, b1 u3 o0 k" z: L
curve_val(2,:)=y_uni;: S0 |. h$ S* T& b/ t
%% //对提取出的数据进行拟合(按实际情况进行修改)$ Z: Q  X$ Z$ r" ^# v
[p,s]=polyfit(curve_val(1,:),curve_val(2,:),4);  %//多项式拟合(为避免龙格库塔,多项式拟合阶数不宜太高)
5 i% i5 z7 X5 |; I6 h1 T) @[y_fit,DELTA]=polyval(p,x_uni,s);  %//求拟合后多项式在x_uni对应的y_fit值2 Y- Z( P' A) P3 z) k
figure,plot(x_uni,y_fit),title('拟合后的曲线')& S3 B4 s) g' o7 e2 W. b
axis([min_x,max_x,min_y,max_y])  %//根据输入设置坐标范围
( |% C) F& B3 ~/ r% {/ k
作者: SsaaM7    时间: 2021-1-20 19:01
利用matlab从图片中提取曲线坐标数据




欢迎光临 EDA365电子论坛网 (https://bbs.eda365.com/) Powered by Discuz! X3.2