帧间预测运动估计：块匹配、光流与凸优化MATLAB实现

简介：这份凸优化大作业代码包面向正在学习凸优化、计算机视觉与视频处理的高校学生，围绕运动估计与帧间预测这一典型课题展开，重点解决如何用多种算法完成 mv 子视角预测并做横向对比的问题。包内共 300 个文件，包含 60 个 m 脚本、4 个 py 脚本、8 个 mat 数据文件以及 216 张 png 结果图，另有若干 txt 说明与 asv 备份，压缩包约 214.99MB，其中 m 与 py 文件承载算法实现，mat 与 png 则保存序列数据与实验输出。方案覆盖运动估计的穷尽块匹配、三步搜索法与分层块匹配三条路线，并延伸出光流估计与单应性矩阵两大类，共组合出十一种方法，每种方法均有独立可运行的代码，便于读者逐项复现、对比误差与耗时差异。受上传体积限制，数据集仅保留少量序列，但仍足以完成完整实验闭环。已有 301 人学习下载，适合借此理解多种预测思路的实现细节与组合技巧。

1. 帧间预测不是“找最像的块”，而是解一个带正则的优化问题

很多人第一次跑完块匹配都有同一个困惑：运动矢量场像一张噪点图，相邻块的方向毫无关联，PSNR 却只涨了零点几 dB。问题往往不在搜索策略，而在于一开始就把运动估计理解成了“在参考帧里找最像的块”，而不是“解一个带正则项的能量泛函”。凸优化这套语言的价值就在这里——它把块匹配、光流法、单应性矩阵这些看起来八竿子打不着的做法，统一写成同一个泛函在不同形变模型、不同罚项下的实例，让你能说清每一步迭代究竟在让什么量下降。这份 MATLAB 实现的十一种方法，恰好把这根轴完整摊开了一遍。

2. 从块匹配到光流：十一种方法其实在解同一个能量泛函

把十一种方法并排放进一个工程里，最容易产生的错觉是“这是一堆互不相干的技巧集合”。恰恰相反，它们的差别只有三个：形变模型有多自由、数据项用什么范数度量、正则项怎么约束。把这三件事拆开，十一个名字立刻塌缩成一张二维表。

2.1 数据项与正则项：所有方法共用的那一个泛函

几乎所有帧间预测方法，无论包装成什么形式，最终都在最小化这样一个东西：

E(u, v) = ∫Ω (I₁(x + u(x)) − I₀(x))² dx + λ ∫Ω ρ(∇u, ∇v) dx

第一项叫数据项，来自亮度守恒假设——同一个物理点在相邻两帧里的灰度不变。它在真实视频里必然不成立：遮挡会让点消失，镜面反射会改变亮度，自动曝光会让整帧整体偏移。所以数据项永远有残差，而残差最大的地方通常正是运动边界，这也是所有方法都躲不开的误差源。

第二项叫正则项，作用是把“残差大”这件事的锅分摊到空间邻域上，避免每个像素各自为政地乱猜。λ 控制分摊强度：λ 趋近 0，解退化成逐像素的独立最小二乘，噪声被完整放大；λ 趋近无穷，解被强压成常数场，运动信息被抹平。真正的手艺全在中间的取值区间里。

不同方法在这两项上的选择，构成了一张很有信息量的对照表：

这张表最有意思的地方在“未知量规模”那一列：从每块 2 个未知量，一路到每帧 2×H×W 个未知量。这其实是一条连续的自由度谱系，而十一种方法只是在这条谱系上取了不同的采样点。

2.2 参数化与非参数化：十一种方法真正的分界线

参数化方法假设整幅图或者每个块的形变可以用少量参数描述。平移是 2 个参数，仿射是 6 个，单应性是 8 个，双线性 warp 是 8 个。这类方法的好处是约束极强，外点几乎不可能把解带偏——只要参与拟合的点里超过一半是内点，最小二乘解的方向就大致正确，再叠一层 RANSAC 做剔除基本能收敛。坏处是它无法表达非刚体运动：一个人的手臂摆动、一片树叶的抖动，用全局单应性去拟合只会得到一张平滑的、谁都不像的图。

非参数化方法则给每个像素分配独立的位移，未知量规模直接爆炸到百万级。这时候如果没有强正则项，问题就是严重欠定的——方程数（每像素 1 个亮度约束）远小于未知量数（每像素 2 个位移）。光流法能成立，全靠正则项把欠定问题补成适定问题。Horn-Schunck 用 L2 平滑项换来闭式迭代，代价是运动边界被抹平；TV-L1 用全变差换来保边，代价是要解一个非光滑的鞍点问题。

块匹配处在两者之间，很少有人这么描述它。它其实是“分段常数”的参数化：每块内部强制位移恒定，块与块之间完全不耦合。所谓块效应，本质就是形变模型在块边界处不连续，跟压缩带来的块效应是同一个数学来源。网格 warp 则是把这种分段常数换成分段线性，节点之间做插值，于是块边界的不连续变成了连续的折线，这也是它在大形变场景下比块匹配舒服得多的原因。

理解了这条谱系，再看“十一种方法”就不会觉得是罗列。它们从强约束弱表达到弱约束强表达，构成了一个完整的取舍空间，而每一档取舍都对应一类具体场景。

2.3 凸性从哪来：线性化之后的子问题才是凸的

这是整套东西里最容易被跳过、但最该讲清楚的一点。

原始的数据项 (I₁(x + u) − I₀(x))² 关于 u 是高度非线性的——你把图像当成 u 的函数，它是一堆局部极值堆出来的曲面，凸性根本无从谈起。任何声称“用凸优化解运动估计”的说法，都需要补一个前提：先做一阶线性化。

把 I₁(x + u) 在 u = 0 处展开，取一阶项：

I₁(x + u) ≈ I₁(x) + Ix·u + Iy·v

其中 Ix、Iy 是参考帧的空间梯度，把 I₁(x) − I₀(x) 记作 It。数据项就变成了

(Ix·u + Iy·v + It)²

这是关于 (u, v) 的二次型，半正定，天然凸。整个问题凸不凸，此时只取决于正则项 ρ 是不是凸函数：L2 是凸的，L1 是凸但不可微，全变差是凸但更不可微，而截断二次、Lorentzian 这类鲁棒罚项是非凸的，一旦用上就得靠初值或者由粗到细的延拓策略来绕。

所以准确的说法是：凸优化在这里提供的不是“全局凸的目标函数”，而是“每次线性化后都有一个凸的子问题可以可靠求解”。这也解释了为什么金字塔是必需品而不是优化项——线性化只在位移远小于一个像素时成立，位移一大，一阶近似就彻底失效，子问题再凸也解不到正确的方向上去。

3. 运动估计的实现骨架：MATLAB 里那些矩阵到底在算什么

理论写到这儿，剩下的问题全是“矩阵里的元素到底来自哪”。这一章把三种代表性方法的核心循环摊开看。

3.1 块匹配：代价函数与搜索策略的代码落点

全搜索是所有方法里最好写、也最好拿来当基准的。它的核心就两个循环，但每一行都藏着口径问题：

function [mv, cost] = fullSearchBlock(ref, cur, bx, by, B, S) % ref/cur : 单通道 double 图像，已归一化到 [0,1] % (bx,by) : 当前块左上角坐标，注意是 (列, 行) % B : 块边长；S : 搜索半径 patch = cur(by:by+B-1, bx:bx+B-1); % 当前块，作为模板 best = inf; mv = [0 0]; for dy = -S:S for dx = -S:S cand = ref(by+dy:by+dy+B-1, bx+dx:bx+dx+B-1); % SAD 代价，等价于 L1 数据项；换成 .^2 求和就是 SSD(L2) d = sum(abs(cand(:) - patch(:))); if d < best best = d; mv = [dx dy]; % 记录的是整数像素位移 end end end cost = best;

bx:bx+B-1这种切片写法隐含了一个约定：块坐标和图像索引用的是同一套原点。一旦你在别处把坐标定义成像素中心、这里却用整数索引，半像素的偏移就会悄悄出现，后面插值全歪。

代价函数的选择比想象中重要。SAD 对应 L1 数据项，对光照的整体偏移不敏感，且没有乘法，硬件上极友好；SSD 对应 L2，对高斯的传感器噪声是最大似然解，但对个别异常像素极度敏感——一个坏点就能把整块的匹配结果拉走。这也是为什么实际编码器普遍用 SAD 加一个小的偏置项，而不是直接上 SSD。

搜索策略那一层更值得说清：三步搜索、菱形搜索、六边形搜索，它们不是优化算法，而是启发式采样。整数格点上的代价曲面是一个布满局部极小的离散曲面，梯度信息完全不可用，能做的只有“少采几个点、赌最优点在上一轮最优点附近”。赌对了省 90% 的计算量，赌错了就掉进局部极小，PSNR 掉零点几 dB 悄无声息。

整数精度定下来之后还有亚像素。常见做法是在最优点附近对代价做一维抛物线拟合，或者直接对参考帧做插值滤波后重跑一遍。前者便宜但精度有限，后者贵一些但和编码器里的插值滤波器口径一致，后续做对比实验时不会出现“我明明记得用了同样的插值”这种对不上的情况。

3.2 光流：Horn-Schunck 与 Lucas-Kanade 差在哪个矩阵上

把 2.3 节那个二次型加上 L2 平滑项，写出欧拉-拉格朗日方程，离散化后得到的迭代式就是 Horn-Schunck：

function [u, v] = hornSchunk(I0, I1, alpha, iters) [Ix, Iy] = gradient(I0); % 中心差分，注意 gradient 返回的是 d/dx, d/dy It = I1 - I0; % 时间梯度：后帧减前帧 u = zeros(size(I0)); v = u; K = [1 1 1; 1 0 1; 1 1 1] / 12; % 邻域均值核，中心为 0——离散拉普拉斯的关键 for k = 1:iters ubar = conv2(u, K, 'same'); % 邻居的平均，不含自己 vbar = conv2(v, K, 'same'); denom = alpha^2 + Ix.^2 + Iy.^2; % 每个像素解一个 2x2 正规方程，闭式写出 u = ubar - Ix .* (Ix.*ubar + Iy.*vbar + It) ./ denom; v = vbar - Iy .* (Ix.*ubar + Iy.*vbar + It) ./ denom; end

K这个核是全篇最容易写错的地方。离散拉普拉斯在欧拉-拉格朗日方程里贡献的是 (邻居均值 − 自身)，所以核的中心必须是 0，周围六个数加起来除以 12。有人图省事写成ones(3)/9，中心不为零，迭代就变成了在解另一个方程，收敛结果整体偏软，运动边界糊成一片，而且怎么调 alpha 都救不回来。

参数含义也值得逐个数一遍。alpha是平滑强度的平方根量级，它和梯度幅值的量纲必须匹配：图像归一化到 [0,1] 时，梯度典型量级在 0.01 到 0.1，alpha取 0.5 到 2 都合理；如果图像还是 0-255 的 uint8 直接转 double 没归一化，梯度能到几十，alpha不放大到几十倍，平滑项就等于不存在，迭代几步就发散。iters不是越大越好，Horn-Schunck 的迭代等价于对稀疏线性系统做 Jacobi/Gauss-Seidel，收敛速度是 O(1/k)，几百次之后改进微乎其微，而误差早就传播到整幅图了。

Lucas-Kanade 走的是另一条路：不假设全局平滑，而是假设一个小窗口内位移恒定，然后对窗口内所有像素的方程做最小二乘。窗口一求和，得到的就是结构张量：

Sxx = conv2(Ix.^2, k, 'same'); Sxy = conv2(Ix.*Iy, k, 'same'); Syy = conv2(Iy.^2, k, 'same'); Sxt = conv2(Ix.*It, k, 'same'); Syt = conv2(Iy.*It, k, 'same'); detS = Sxx.*Syy - Sxy.^2; u = (-Syy.*Sxt + Sxy.*Syt) ./ detS; v = ( Sxy.*Sxt - Sxx.*Syt) ./ detS;

真正的分水岭在detS。当窗口内的梯度全都指向同一个方向——比如一条纯水平的边缘——结构张量就退化成秩 1 矩阵，detS趋近 0，除法把数值噪声放大成天文数字。这就是孔径问题的矩阵形式。工程上的标准做法不是加一个 epsilon 硬顶，而是算结构张量的两个特征值，用最小特征值做门限：λmin 小的位置直接判为不可靠，输出零位移或者交给插值补全。这一步做不做，稀疏光流的轨迹质量差得非常明显。

3.3 单应性矩阵：从 DLT 到 SVD 零空间，再到逆向 warp

单应性是 8 自由度的射影变换，四对不共线的点就能唯一确定。写成齐次形式就是一组线性方程 Ah = 0，h 是 9 维向量：

function H = dltHomography(p1, p2) % p1,p2 : N×2，N>=4；p1 是参考帧坐标，p2 是当前帧坐标 [T1, n1] = normalizePoints(p1); % Hartley 归一化 [T2, n2] = normalizePoints(p2); A = zeros(2*size(n1,1), 9); for i = 1:size(n1,1) x = n1(i,1); y = n1(i,2); X = n2(i,1); Y = n2(i,2); A(2*i-1,:) = [-x -y -1 0 0 0 X*x X*y X]; A(2*i, :) = [ 0 0 0 -x -y -1 Y*x Y*y Y]; end [~, ~, V] = svd(A); % Ah=0 的最小二乘解 = 最小奇异值的右奇异向量 Hn = reshape(V(:,end), 3, 3)'; H = T2 \ Hn * T1; % 反归一化回原坐标系 H = H / H(3,3); % 固定尺度，让 H(3,3)=1

每一对点贡献两行方程，N 对点贡献 2N 行、9 个未知量。齐次方程组的零解没有意义，所以约束 ||h||=1，问题变成在单位球面上最小化 ||Ah||，解就是 A 的最小奇异值对应的右奇异向量。这一步的几何含义是：h 落在 A 的近似零空间里，DLT 本质是在找一个让所有点对残差平方和最小的投影矩阵。

normalizePoints那两行看起来多余，实际是整个函数能不能用的关键。原始像素坐标动辄成百上千，A 的列量级相差三个数量级，条件数轻松上 1e10，SVD 返回的“最小奇异向量”基本是数值噪声。Hartley 归一化的做法很朴素：把点平移到质心、再缩放到平均距离为 √2，A 的条件数立刻回到个位数。少了这两行，warp 出来的图像会扭成一团谁也认不出的东西，而且报错信息极其含糊。

拿到 H 之后做重采样，必须用逆向映射：

function out = warpInverse(img, H, outSize) % 对输出图像的每个像素求它在输入图中的来源位置，避免出现空洞 [xx, yy] = meshgrid(1:outSize(2), 1:outSize(1)); p = H \ [xx(:)'; yy(:)'; ones(1, numel(xx))]; % H 把输入坐标映到输出坐标 src = p(1:2,:) ./ p(3,:); % 齐次坐标归一化 out = interp2(img, src(1,:), src(2,:), 'linear', 0); out = reshape(out, outSize);

正向映射会把多个输入像素打到同一个输出位置，剩下的位置没人填，图上出现密密麻麻的麻点；逆向映射保证输出每个像素都被访问且只被访问一次，interp2的双线性插值顺带解决了亚像素取值。最后一个参数0是越界填充值，用它还是用边缘复制，会直接影响后续计算 PSNR 时的边界区域——这是一个纯粹的口径问题，但足以让两组实验对不上。

4. 求解器层：怎么把预测误差写成标准凸形式

前两章解决的是“问题长什么样”，这一章解决“交给谁去解”。同一个能量泛函，罚项换一下，可用的求解器和最终画质都会变。

4.1 罚项选择：L2、L1、TV 分别把问题交给谁

线性化之后的数据项是二次型，可以写成 ||Au − b||²，A 的每一行是一个像素的 [Ix, Iy]。剩下的就是正则项怎么选。

L2 罚项让整个问题变成一个严格凸的二次规划，正规方程是 (AᵀA + λL)u = Aᵀb，其中 L 是离散拉普拉斯矩阵。这个系统是对称正定的稀疏系统，共轭梯度几十次迭代就能收敛到很高的精度，而且解唯一、对初值不敏感。代价是 L2 平方惩罚大残差，导致运动边界——那里残差本来就是大的——被强行压低，结果是整幅光流场又糊又软，边界向内收缩。Horn-Schunck 的所有优点和所有缺点都来自这一个选择。

L1 罚项换成绝对值，问题仍然凸但不可微。常见解法是 IRLS，把 |x| 近似成 x²/|x| 迭代加权，内层还是解二次规划；或者直接用 ADMM 拆成“二次项 + 软阈值”两步。L1 的好处是允许少量大残差存在，遮挡区域和运动边界不会被过度平滑，代价是解可能不唯一，迭代次数明显变多。

全变差把罚项施加在梯度的模上，既保边又不像各向异性 L1 那样产生阶梯状伪影。代价是要解一个鞍点问题，Chambolle-Pock 那套原始对偶迭代是标准解法：交替做一次对偶上升和一次原始下降，每步都是逐像素的闭式更新，实现起来意外地短，但步长和对偶变量的初值需要调，步长取大了直接发散。

选择哪种罚项，本质上是在问“我愿意为保边付多少计算代价”。对帧间预测这个任务来说，多数情况下 L2 已经够用，因为最终目的不是拿到漂亮的光流场，而是让重建帧的残差最小；过度追求光流场的视觉质量，反而可能在 PSNR 上吃亏。这是一个很容易被忽略的取舍：优化目标和评价指标不是同一个东西。

4.2 拉格朗日乘子与率失真：λ 到底是什么东西

视频编码里的帧间预测还有一个额外约束——运动矢量本身要编码进码流，占用比特。这就把问题从单纯的“最小化预测误差”变成了“最小化误差与比特的加权和”：

J = D + λ·R

D 是预测残差的失真（常用 SAD 或 SSD），R 是编码这个运动矢量所需的比特数估计。这个式子的形式眼熟，因为它就是一个拉格朗日函数：在给定比特预算下最小化失真，等价于对某个 λ 最小化 J。

λ 的几何含义是率失真曲线上的负斜率。把每种编码方案画成 (R, D) 平面上的一个点，所有可行方案的下凸包就是理论最优边界，λ 就是这条边界在某点的切线斜率。给定一个 λ，最小化 J 得到的点必然落在下凸包上；反过来，凸包上方的任何点，无论怎么调参数都取不到——因为它是被凸包支配的次优解。工程上常见的做法是让 λ 跟量化步长成幂律关系，量化越粗，λ 越大，编码器越倾向于用更少比特的粗糙运动。

这条结论直接解释了为什么“十一种方法”放在一起有价值。不同的运动估计方法在不同的比特区间上才有优势：块匹配在低比特、大位移场景下无可替代，因为它的运动矢量编码代价低且搜索范围可控；稠密光流在高比特、精细形变场景下才体现出补偿精度优势，但把逐像素的位移编进码流的代价高得离谱，必须靠参数化模型压缩。如果只实现一种方法，你看到的只是率失真曲线上的一个点；把十一种都跑一遍，才看得到那条凸包的形状，也才知道自己手里这个方法到底卡在哪个位置。

4.3 由粗到细：金字塔救的不是凸性，是线性化

金字塔常被当成“加速技巧”，这个理解偏差很大。它的真正作用是扩大线性化的有效范围。

一阶泰勒展开成立的隐含条件是位移在一个像素量级以内。位移是 10 个像素时，Ix·u 这一项已经完全不能代表真实的灰度变化，线性化后的二次型指向的方向和真实下降方向可能相差几十度，子问题再凸也只会稳稳地收敛到错误的位置。

金字塔的做法是把图像逐级降采样，位移也随之按比例缩小。10 个像素的位移在四分之一分辨率上只剩 2.5 个像素，在十六分之一分辨率上不到 0.7 个像素，线性化重新变得可信。然后在粗层解出一个粗略的位移场，上采样作为细层的初值，细层只需估计残差位移。

这在优化理论里叫延拓法或者连续化方法：把一个难解的非凸问题，换成一串好解的近似问题。每一次上采样都是一次热启动，让细层的迭代从一个已经很接近真值的位置出发。这解释了实践中的几个现象：不做金字塔的 LK 在大位移上完全失效，做三层金字塔就能追上；金字塔层数不是越多越好，层数过多时最粗层本身已经看不清纹理，估出来的初值反而把细层带偏；以及为什么粗层可以多迭代几次、细层少迭代——粗层的解决定了大方向，细层只是打磨。

5. 常见问题排查：那些让 PSNR 莫名掉 3dB 的坑

这套东西跑通不难，跑对很难。下面几条都是“代码没报错、图也能出、结果就是不对”的典型情况。

5.1 坐标与插值口径

现象：运动矢量场全是 (0, 0) 或者全都贴在搜索窗的边界上，重建帧的 PSNR 比直接不做预测还低一点点。

原因：坐标约定混用。MATLAB 的矩阵索引是先行后列，而图像坐标习惯写成 (x, y)，两者一旦在某个函数里互换，块的切片位置就整体转了 90 度，匹配到的自然是完全不相关的区域。另一种更隐蔽的情况是像素中心偏移：整数索引 i 对应的物理坐标是 i，但像素中心在 i+0.5 附近，做亚像素插值时不减这 0.5，误差恒定为半个像素，表现为整体平移方向的系统性偏差。

解决：全工程统一用 (列, 行) 的切片写法，只在几何变换的入口和出口做一次显式的坐标转换并写注释。做亚像素插值前先把坐标减 0.5，插完再加回来。做完之后拿一张纯平移的合成图验证，估计出来的位移应该精确等于真值，误差在 1e-6 量级，多一位小数都说明还有口径问题。

5.2 数值与求解器

现象：光流迭代到第十几次，u、v 的数值开始指数上涨，几十次后变成 NaN，图上出现大片纯白。

原因：图像没有归一化。uint8 转 double 之后取值范围仍然是 0 到 255，梯度幅值能到几十，数据项比平滑项大四个数量级，denom里的 alpha² 完全不起作用，迭代等价于在解一个病态系统。另一个常见诱因是帧差 It 的符号反了，It 取反等于把下降方向整体调头，梯度法在这种情况下不是不收敛，而是朝反方向稳步发散。

解决：所有进入求解器的图像先除以 255 归一到 [0,1]，alpha 按梯度能量的均方根做缩放而不是拍脑袋定值。迭代过程中加一个断言，检查 max(|u|) 是否超过搜索范围的合理上限，超了立刻报错而不是让它继续算——早失败比晚失败省时间。

5.3 数据与评价口径

现象：和参考实现对比，PSNR 总是稳定差 1 到 2 dB，逐帧检查也找不到明显异常。

原因：三个地方最常出问题。一是颜色空间：把 RGB 三个通道分别估计光流再拼回去，三个通道的梯度方向不一致，会得到三个不同的位移场，合成后边缘出现彩色镶边；正确做法是转到亮度色度空间，只在前者上估一次，三通道共用。二是边界处理：conv2的'same'参数默认补零，图像四周一圈的梯度是被零值污染过的，边界区域的光流不可信，如果评价时把这一圈也算进去，分数就被拉低。三是参考帧来源：编码器里重建帧必须用解码后的重建帧，用原始帧当参考会得到虚高的 PSNR，一对接真实码流就露馅。

解决：固定一条处理链并把中间量全部落盘——归一化后的图像、梯度图、代价曲面、光流场、warp 结果各存一份。对比时逐级看，第一个出现差异的环节就是问题所在。这套“中间量全留”的习惯非常土，但每次都能把定位时间从半天压到十分钟。

6. 设计取舍：十一种方法放在一起，价值不在“挑一个最好的”

把十一种方法塞进同一个工程，最没用的读法就是跑一遍 benchmark、挑分数最高的那个用。它们真正的价值在于共用了一套数据接口和评价口径，从而让消融变得可信。

取舍的核心是三个量的拉锯：形变模型的自由度、正则项的强度、求解器的计算量。自由度上去，能表达的运动更复杂，但需要更强的正则来压制欠定，计算量也跟着涨；正则加强，数值更稳，但边界和细节被抹平；求解器换成保边的罚项，画质好看，迭代次数翻几倍。任何一次调参都在这三者之间挪动，不存在同时最优的点。

验证方法我一般这么设计：先用合成的已知变换做回归测试——给一张图施加已知的平移、旋转、单应性变换，看估计出来的参数和真值的误差是否在数值精度内。这一步能把绝大多数坐标和插值口径问题挡在门外。再用真实序列看率失真曲线，把十一种方法的 (R, D) 点画在同一张图上，看谁落在凸包上、谁被支配。落不到凸包上的方法不是“实现错了”，而是它的适用区间不在这里——这个判断比单纯比较 PSNR 有用得多。

有个坑我踩过不止一次：调完参数发现某一帧的 PSNR 突增，兴奋半天，最后发现是那一帧恰好是 I 帧，根本没走帧间预测路径。从那以后我每次做方法对比，都强制先把帧类型、参考帧索引、量化参数一起打印出来对齐，确认对比的是同一件事再谈数字。希望帮到你。

本文还有配套的精品资源，点击获取