★ FEATURED ARTICLE 简介面向科研人员与MATLAB学习者的曲线图片识别工具包可从论文截图或实验图像中提取曲线坐标免去重复编写绘图代码便于多来源数据对比与复核。压缩包共14个文件、约459KB以12张PNG示例图片、1个MAT文件与1个M脚本组成图像样本用于测试识别流程MAT数据与脚本便于加载结果、复现实验并继续调整算法。程序涉及灰度化与二值化、Canny等边缘检测、连通域标记与轮廓追踪、曲线拟合及绘图输出也包含图像读写和结果保存等基础操作适合具备一定MATLAB基础、希望把图像处理与数据分析串起来的中级读者。目前已有2019人学习可作为课程实验、论文复现或科研数据整理时的参考实现。1. 从一张截图里把曲线捞回来这个 matlab 程序解决的是“只有图没有数”的麻烦手里只剩一张曲线截图要复现实验就得把数据点一个个描出来——这事我干过不止一次。某开发者拿一份报告里的十几张 matlab 出图来找我原始 .mat 早不知道丢哪儿去了只剩 PNG。这个曲线图片识别程序干的就是这活读入一张曲线图先摘掉网格和边框、锁定绘图区再按你给的标定点把像素坐标换算成数据坐标最后吐出 CSV。它适合两类人一类是手里有图无数据、要还原趋势去做对比的另一类是想拿 matlab 把图像处理流水线走一遍的。前提是图得干净——白底、轴清楚、曲线颜色能和背景分开那种糊成一片的论文扫描件别指望它能自己变魔术。2. 项目骨架与技术选型五个函数串起一条“像素到数据”的流水线这类工具最怕写成一个上千行的巨型脚本改一个阈值要翻半天。我的拆法是按处理阶段切函数每个函数只干一件事、只吃一个 cfg 结构体中间结果老老实实返回出问题能单步停下来看。matlab 在这一点上其实比很多人想的顺手Image Processing Toolbox 里的imbinarize、bwmorph、regionprops、bwlabel已经把二值图那套基本操作封得够用矩阵索引又能让你一列一列扫像素不用像在别的语言里那样先想尽办法把图转成数组再转回来。2.1 目录结构与各模块职责整个工程保持扁平一个入口、一个配置、五个阶段函数外加一个 utils 目录放不进入主流程的小工具。文件职责输入 → 输出main.m唯一入口串流程、抓异常、写日志图片路径 → 无config.m所有阈值、标定值、开关集中在这里无 → cfg 结构体preprocess.m读图、灰度化、对比度拉伸、二值化路径 cfg → bw / gray / rawremove_grid.m形态学去网格线、去边框bw cfg → bwlocate_axes.m投影法锁定绘图区上下左右边界bw cfg → box 结构体color_mask.mHSV 阈值分割按颜色拆出曲线raw cfg → 逻辑掩膜extract_curve.m逐列扫描取中值像素→数据坐标mask box T → x / ypost_process.m去重、排序、平滑、野值剔除x / y cfg → x / yexport_data.m写 CSV 并画回原图做叠图核验x / y → 文件utils/split_segments.m把一列里的行号按间距拆段rows 容差 → cell把阈值全部塞进config.m是血泪经验换来的。早期我把二值化阈值直接写死在函数里换个图源就得全局搜索替换改到第三张图的时候人已经麻了。现在换数据集只动配置文件主流程一行不碰。2.2 依赖清单与版本差异必须的是 Image Processing Toolbox其余都是基础功能。有两个函数要留意版本差异imbinarize是较新版本才有的老版本得退回im2bw配graythreshwritematrix也是后来才加的老版本用csvwrite或dlmwrite代替。我在函数里做了个简单的能力判断避免拿老版本跑的时候第一个函数就报未定义。提示如果打开就提示某个函数未定义先别怀疑代码ver看一眼工具箱装没装、版本够不够这一步能省掉大量无效调试。曲线拟合和插值用的polyfit、interp1、movmean、medfilt1都是基础库自带其中medfilt1属于 Signal Processing Toolbox不想多装一个工具箱就用movmedian代替效果差别不大。2.3 一次调用的主流程 main.m入口函数刻意写得很薄全流程读下来不超过二十行出问题时一眼能定位是哪个阶段挂了。function main(imgPath) % 单张曲线图的识别入口 cfg config(); % 所有阈值和标定值都从这里来 try [bw, gray, raw] preprocess(imgPath, cfg); % 1. 灰度 二值化 catch ME fprintf([读图失败] %s - %s\n, imgPath, ME.message); return; end bw remove_grid(bw, cfg); % 2. 去网格线与外框 box locate_axes(bw, cfg); % 3. 锁定绘图区四条边 if isempty(box.left) fprintf([定位失败] %s 找不到坐标轴检查背景是否太脏\n, imgPath); return; end T calibrate(box, cfg.axisVal, cfg); % 4. 建像素-数据的映射 mask color_mask(raw, cfg); % 5. 按颜色剥出曲线 [x, y] extract_curve(mask, box, T, cfg); % 6. 逐列取点 [x, y] post_process(x, y, cfg); % 7. 去重去野值 export_data(x, y, cfg, imgPath); % 8. 导出 叠图核验 end参数上imgPath支持相对路径和绝对路径cfg.axisVal是一个四元向量[x左 x右 y下 y上]填的是坐标轴上对应的真实数值这个值目前靠人工读图填我在第 4 章会讲为什么暂时不打算做全自动刻度 OCRcfg里另有binThresh、gridLenH、gridLenV、minArea、hRange、sMin几个关键项后面逐个交代。流程里有两处提前 return——读图失败和定位失败。别小看这两个判断批量跑几十张图的时候没有它们你会得到一堆莫名其妙的空 CSV还得回头猜是哪张图出的事。3. 预处理与坐标轴定位先把网格、边框、刻度从图里摘干净曲线提取的准确率七成取决于这一步干不干净。剩下的三成里两成是标定一成才是取点算法。很多人上来就写逐列扫描扫出来的点全是网格线的位置画回去一看是一条锯齿状的鬼画符这就是跳过了预处理。3.1 灰度化、对比度拉伸与自适应二值化先压灰度再做对比度拉伸。截图和扫描件的背景很少是纯白往往是 240 左右的浅灰直接二值化会把大片背景也算成前景后面bwareaopen都救不回来。function [bw, gray, raw] preprocess(imgPath, cfg) raw imread(imgPath); if size(raw, 3) 3 gray rgb2gray(raw); % 彩色图压灰度颜色信息留给 color_mask 用 else gray raw; end gray imadjust(gray); % 拉伸对比度把发灰的背景推到接近纯白 if isempty(cfg.binThresh) % 光照不均时用自适应Sensitivity 越小越保守 bw imbinarize(gray, adaptive, ... ForegroundPolarity, dark, Sensitivity, cfg.sensitivity); else bw imbinarize(gray, cfg.binThresh); % 固定阈值适合同批同源图 end bw bwareaopen(bw, cfg.minArea); % 去掉零散噪点 endForegroundPolarity设成dark是关键它告诉函数“我要找的是比背景暗的东西”曲线和轴线都是深色背景是浅色这个极性判断错了会得到一张反色的图后面所有步骤全废。Sensitivity我一般从 0.4 起调调大到 0.6 以上时噪声会明显变多网格线的断点也容易被连起来反而更难处理。minArea默认给 20 像素图片分辨率越高这个值越要往上加1080p 的图给到 50 到 80 比较合适。固定阈值和自适应阈值的选择有个简单判断标准同一批图如果来自同一个渲染器、背景一致用固定阈值binThresh一般 0.85 到 0.92 之间结果最稳定如果是扫描件、拍照件、或者图里有渐变背景老老实实走自适应。3.2 形态学去网格的边界条件网格线在形态学眼里特征极其明显又长又直要么纯水平要么纯垂直。开运算用一个长条结构元就能把它们单独提取出来再从原图里减掉。function bw remove_grid(bw, cfg) % 长条结构元只保留足够长且方向一致的线段 hLine imopen(bw, strel(line, cfg.gridLenH, 0)); % 水平网格 vLine imopen(bw, strel(line, cfg.gridLenV, 90)); % 垂直网格 grid hLine | vLine; if cfg.keepFrame % 保留外框只减掉框内部的网格 inner false(size(bw)); inner(cfg.frameTop:cfg.frameBottom, cfg.frameLeft:cfg.frameRight) true; grid grid inner; end bw bw ~grid; bw bwareaopen(bw, cfg.minArea); endgridLenH和gridLenV的长度要卡在一个区间里比一格网格长、比整条坐标轴短。设短了网格断点减不干净设长了坐标轴本身也被当成网格一起减掉后面locate_axes直接找不到边。我的做法是先跑一次不管网格的版本locate_axes输出四条边的位置量一下绘图区高度取gridLenV ≈ 绘图区高度 × 0.7这个经验系数在多数图上都成立。keepFrame这个开关背后是一个很实际的取舍如果坐标轴本身就是背景色、只有网格那就得靠外框来定边界如果图是标准的 matlab 出图、四周有框线那就必须保留外框。判断方法很简单把hLine和vLine显示出来看一眼就明白了。这里有个必须提前说清的坑如果某条曲线本身接近水平比如一段慢变的温度曲线imopen会把它当成网格线减掉。所以remove_grid只应该作用在灰度二值图bw上用来找轴真正的曲线提取走的是color_mask的颜色通路两条路互不干扰。这个设计是刻意的不是冗余。3.3 投影法锁定绘图区四条边坐标轴是大面积连续深色像素在行方向和列方向的投影里会形成尖峰。这个方法老套但极其可靠几行代码就能定位。function box locate_axes(bw, cfg) colSum sum(bw, 1); % 每列的前景像素数 - 竖线 rowSum sum(bw, 2); % 每行的前景像素数 - 横线 colSum movmean(colSum, 3); % 轻微平滑压掉单像素毛刺 rowSum movmean(rowSum, 3); thC max(colSum) * cfg.axisRatio; % 只保留接近最高峰的列 thR max(rowSum) * cfg.axisRatio; candC find(colSum thC); candR find(rowSum thR); if isempty(candC) || isempty(candR) box struct(left, [], right, [], top, [], bottom, []); return; end box.left min(candC); box.right max(candC); box.top min(candR); box.bottom max(candR); endaxisRatio默认 0.6含义是“至少要有最高峰六成的长度才算轴”。设太低会把密集的网格候选也吃进来导致box被撑大到整张图设太高则在轴被曲线截断的情况下漏检。图里如果同时存在左右两条竖轴双 y 轴图min和max会把整个绘图区框住这恰好是我们想要的——两条轴的像素范围就是数据范围。movmean的窗口给 3 就够了给大了会把两条靠得很近的轴合并成一条边界直接偏移几个像素标定误差就是这么来的。3.4 刻度线识别与标定点取值到这一步绘图区的像素范围有了但“左上角对应 (0, 100)”这件事程序并不知道必须由外部告诉它。全自动刻度识别听着很美实际做起来要先切出轴标签区域、再做字符分割、再识别数字、还得处理和科学计数法、负号、小数点投入产出比极低而且识别错一个数字整条曲线就全偏了人还很难发现。我采用的是半自动程序在图上画出检测到的四条边人工确认后用ginput点两个点或者直接在cfg.axisVal里填数值。function T calibrate(box, axisVal, cfg) % axisVal [x左 x右 y下 y上]人工读图填入 px [box.left, box.right]; py [box.bottom, box.top]; % 注意像素行号向下增大y 轴向上增大 if strcmpi(cfg.xScale, log) T.px2x polyfit(log10(px), log10([axisVal(1) axisVal(2)]), 1); T.xMode log; else T.px2x polyfit(px, [axisVal(1) axisVal(2)], 1); T.xMode linear; end if strcmpi(cfg.yScale, log) T.py2y polyfit(log10(py), log10([axisVal(3) axisVal(4)]), 1); T.yMode log; else T.py2y polyfit(py, [axisVal(3) axisVal(4)], 1); T.yMode linear; end T.box box; T.axisVal axisVal; end用polyfit做一次多项式拟合本质就是两点连线但写成拟合的形式后面好扩展——如果你能提供三个以上的已知点比如网格交点直接换成最小二乘抗标定误差的能力立刻上一个台阶。axisVal允许左大右小也就是反向轴polyfit会自动得到负斜率不用特殊处理这点比手写线性插值省心。对数轴的判断必须由使用者显式告知程序没法从像素里看出这是对数轴。判断依据其实肉眼就能看出来刻度是 1、10、100、1000 这种等比排列的就是对数轴。4. 标定与曲线提取把像素坐标反解成物理量坐标轴定好了剩下的核心问题只有一个给定一个像素点 (row, col)怎么算出它对应的 (x, y)。这一步没有玄学就是两次一维映射但坑都藏在细节里。4.1 线性标定与对数轴、反向轴的处理标定公式在第 3 章已经给出来了这里补一句取值时的写法因为对数轴和非对数轴的取反操作不一样混用会得到一列负数或者 NaN。function v px2val(p, T, dim) % dim x 或 y统一入口避免到处写 if if strcmpi(dim, x) f T.px2x; mode T.xMode; else f T.py2y; mode T.yMode; end if strcmpi(mode, log) v 10 .^ polyval(f, log10(p)); % 先取对数域的值再还原 else v polyval(f, p); end endpolyval比polyfit更适合在循环里反复调用它是纯求值没有拟合开销。对数轴这条分支里log10(p)要求像素坐标为正所以box的边界值不能出现 0出现 0 说明定位结果本身是错的先去查locate_axes。反向轴的情况上面已经覆盖了axisVal(1) axisVal(2)时拟合出来斜率为负取值自然就是反向的不需要额外加if。4.2 颜色分割剥出多条曲线一张图里三条曲线是常态靠灰度值区分基本不可能——它们的灰度往往一样深只有色相不同。所以曲线提取这条路走的是 HSV 空间。function mask color_mask(raw, cfg) if size(raw, 3) 3 % 灰度图没有颜色可用退回亮度阈值效果会打折 gray im2double(raw); mask gray cfg.grayCurveThresh; return; end hsvImg rgb2hsv(raw); h hsvImg(:,:,1); s hsvImg(:,:,2); v hsvImg(:,:,3); if cfg.isRedHue % 红色在色相环上跨 0 边界必须用两段拼 hueOK (h cfg.hRange(2)) | (h cfg.hRange(1)); else hueOK (h cfg.hRange(1)) (h cfg.hRange(2)); end mask hueOK (s cfg.sMin) (v cfg.vMin); mask bwareaopen(mask, cfg.minArea); mask imclose(mask, strel(disk, cfg.closeR)); % 补上抗锯齿留下的断点 endsMin是饱和度下限我一般给 0.3 到 0.5。抗锯齿的边缘像素饱和度会明显下降sMin设太高会让曲线变细甚至断开设太低则会把浅色的网格和背景杂色一起吃进来。vMin是亮度下限用来排除黑色文字和坐标轴——文字和轴通常是低饱和度的深色饱和度这一关基本就能挡住vMin是双保险。closeR做一次闭运算半径给 1 到 2 像素专门补抗锯齿造成的单像素缺口。给大了曲线会明显变粗还会把靠得近的两条曲线粘成一坨这个半径一定要小。4.3 列扫描取中值解决线宽、锯齿与自交有了掩膜逐个列扫过去就行。每一列可能命中一片连续的像素行也可能命中好几段——曲线自交、或者同色的两条线交叉时就会出现多段。取中值能抵消线宽带来的偏移取均值则会被抗锯齿的灰边拉偏。function [x, y] extract_curve(mask, box, T, cfg) x []; y []; for c (box.left cfg.margin) : (box.right - cfg.margin) rows find(mask(:, c)); if isempty(rows), continue; end segs split_segments(rows, cfg.gapTol); % 一列内按间距拆段 for k 1:numel(segs) r median(segs{k}); % 中值抵消线宽 x(end1, 1) px2val(c, T, x); y(end1, 1) px2val(r, T, y); end end end function segs split_segments(rows, gapTol) segs {}; if isempty(rows), return; end cut find(diff(rows) gapTol); % gapTol 一般给 2~3 像素 s [1; cut 1]; e [cut; numel(rows)]; for i 1:numel(s) segs{i} rows(s(i):e(i)); end endmargin的作用是把贴着边框的那几列跳过边框本身是深色的很容易被颜色掩膜误吃进来跳过 2 到 3 列就能避开。gapTol控制“多近算同一段”给 2 到 3 像素足够跨过抗锯齿的小缺口给到 8 以上就会把靠得近的两条曲线当成一条。这段代码对一条曲线会输出一组点对 k 条曲线会输出 k 组点但混在一个数组里。实际使用时我会在extract_curve外面套一层按segs的段数做过滤如果某列出现两段以上说明这里有自交或者多条曲线要么按颜色再细分要么在post_process里做轨迹跟踪。4.4 断点、野值与平滑的取舍原始点里一定混着噪声图例框的边缘、曲线穿过网格线的位置、抗锯齿造成的偏移。处理顺序有讲究先排序去重再剔野值最后平滑。function [x, y] post_process(x, y, cfg) [x, idx] sort(x); y y(idx); % 同一个 x 上出现多个点多值曲线时只保留最靠上的那个 [u, ia] unique(round(x, cfg.xDig), last); x u; y y(ia); % 中值滤波剔野值窗口必须奇数 y medfilt1(y, cfg.medWin, omitnan); % 移动均值平滑窗口不要超过曲线最短周期的 1/10 if cfg.smoothWin 1 y movmean(y, cfg.smoothWin); end endxDig是 x 保留的小数位数因为浮点映射出来的 x 几乎不会重复不取整的话unique一个都去不掉。medWin建议 5 或 7专门对付那种突然跳出去一大截的单点野值——它们通常来自图例边框被误识。smoothWin我用 5 到 9再大就会把真实的小幅振荡一起抹平画回原图做叠图对比时一眼就能看出来。平滑是最后一刀做完之后一定要把结果叠回原图上看别只盯着 CSV 里那串数字。叠图核验在export_data里做把提取出的点用另一种颜色画在raw上点是落在曲线上还是漂在曲线旁边肉眼一秒钟就能判断比看残差图快得多。5. 避坑与排查这几处不处理好导出的数据基本不能用下面四条是我在不同图上反复翻车之后记下来的每条都按现象、原因、解决三段走。遇到提取结果不对先按这四条对一遍八成能定位。5.1 曲线整体“变胖”导出值恒定偏移现象提取出的 y 值比真实值稳定偏大或偏小一个小常数曲线形状完全对就是整体上下平移了一点点。原因抗锯齿。渲染器画线时会在主线条两侧生成半透明的过渡像素颜色掩膜的阈值一松这些过渡像素就被算进曲线median取的是加粗后那一坨的中值自然偏移了半个到两三个像素。解决先量线宽再决定阈值。用一行代码看当前掩膜的平均粗细stats regionprops(mask, MinorAxisLength); fprintf(平均线宽 %.2f 像素\n, median([stats.MinorAxisLength]));正常渲染的曲线线宽在 1.5 到 3 像素之间如果量出来 5 以上说明阈值太松。把sMin往上加 0.1、vMin往上加 0.05再量一次。另一种更稳的做法是做骨架化bwskel(mask)把粗线细化成单像素中心线再逐列取点偏移直接消失代价是交叉处骨架会分叉要配合gapTol一起调。5.2 图例和标注被当成曲线现象在 x 轴的某个区间里突然出现一小撮离群点画回原图一看位置正好对图例里那条彩色小线段或某个文字标注。原因图例里的示例线条和曲线同色同宽颜色掩膜分辨不出来文字标注如果是彩色的同理。解决最省事的办法是在config.m里加一个排除矩形cfg.excludeRect [x1 y1 w h]用mask(y1:y1h, x1:x1w) false挖掉。矩形坐标不用手填在图上ginput框两下就出来了。如果想自动一点可以统计每个连通域的包围盒长宽比图例线段通常是短而粗的一小段长宽比接近 3 到 8而真实曲线横跨大半个绘图区长宽比往往在 20 以上用regionprops的BoundingBox一筛就能滤掉大部分。文字标注则靠面积和长宽比双重条件字符的包围盒接近正方形很好认。5.3 对数轴用了线性标定现象曲线形状被压成一条近乎水平的直线或者低值区间的点全挤在一起高值区间拉得特别开。原因图是对数轴cfg.xScale或cfg.yScale忘了设成logcalibrate走了线性分支映射关系完全错位。解决看图判断刻度标签是 1、10、100、1000 这种等比排列的就是对数轴是 0、50、100、150 的就是线性。calibrate里已经按cfg分支处理改配置就行。改完一定要叠图核验因为对数轴的错误映射看起来“也挺像一条曲线”光看 CSV 不容易发现。5.4 曲线自交处一列出现两段现象提取结果在某个 x 附近出现两个 y且两点距离很近画出来是一个小小的竖直跳变。原因曲线自交或者两条同色曲线交叉同一个像素列上确实对应两个数据点。split_segments把它们拆成了两段两段都被写进了结果数组。解决看你的需求。如果只是要趋势post_process里那条unique(..., last)会保留其中一个简单粗暴但够用。如果要还原完整轨迹就得做跟踪以上一列的点为基准在两段里选距离最近的那段继续另一段另起一条轨迹。% 简易跟踪候选段里选离上一列最近的那段 if ~isempty(lastRow) [~, best] min(abs(cellfun(median, segs) - lastRow)); segs segs(best); end这段逻辑放进extract_curve的列循环里就能跑把lastRow在每列末尾更新成当前选中的中值即可。代价是遇到真正的分叉两条曲线从同一点分开会跟丢一条那种情况只能靠颜色或者人工干预别指望纯几何方法能自动分开。6. 二次开发从单张手调改成批量流水线单张图跑通之后真正的需求几乎一定是“我这儿还有两百张”。改造的核心思路是把main.m里那些硬编码的路径和人工交互全部外置让主流程变成一个纯函数输入路径、输出结果中间不弹窗、不停顿。第一步是把配置对象化。config.m改成接受一个名字参数按图源返回不同的 cfgcfg config(datasetA)每批图的阈值、标定值、排除矩形都写在自己的分支里。同批同源的图标定值往往是完全一样的只要绘图区的像素范围一致就能复用T省掉大量重复劳动。第二步是批量入口。用dir收集图片parfor并发前提是每张图互不依赖且不开ginput失败的不中断整个批次单独记到一个失败清单里。files dir(fullfile(cfg.batchDir, *.png)); failList {}; parfor i 1:numel(files) try main(fullfile(files(i).folder, files(i).name)); catch ME failList{end1} sprintf(%s | %s, files(i).name, ME.message); end end第三步是加质检。批量的敌人不是报错是“跑完了但结果不对还没人发现”。我加了一个简单的质量分提取点数和绘图区宽度之比低于 0.3 的判为可疑曲线断得太厉害最高点和最低点落在绘图区边界 2 像素以内的判为可疑很可能是轴被当成了曲线。可疑的图自动进复核队列人工只看这一小撮效率比全量看高一个数量级。检查项阈值判为可疑的含义有效列占比 0.3曲线断裂严重或颜色阈值过严端点贴边距边界 2 px坐标轴被误识别为曲线点数突变与上批均值偏差 50%图源变了阈值需要重调最后一步是把拟合接上去。原始点导出 CSV 之后如果要的是趋势方程而不是逐点数据可以在export_data里追加一段polyfit把阶数和 R² 一起写进结果表方便一眼判断拟合质量。从那以后我每次跑批量识别都强制先把前十张图的标定结果和叠图核验挨个看一遍确认没有系统性偏移再放开全量跑——这个习惯让我躲过了至少两次“两百张图全偏了半个像素只能从头再来”的返工。希望帮到你。本文还有配套的精品资源点击获取
阅读完成 · 觉得有帮助？
