十年匠心定制 · 商业建站与技术教学双线并行 咨询热线:400-886-1026 service@lmnt.cn
ARTICLE DETAIL

资讯详情

深耕网站建设与运营推广的一线实战洞察。

用Python实现EVM算法:放大视频中看不见的微小运动

用Python实现EVM算法:放大视频中看不见的微小运动 如果你盯着一台固定机位的监控画面并且知道画面里某个对象的移动幅度只有一两个像素绝大多数人都会觉得“这跟静止没区别”。EVM算法Eulerian Video Magnification欧拉视频放大就是为了对付这类问题的典型方法不靠光流跟踪不靠特征点匹配而是把视频里肉眼根本察觉不到的微小运动和颜色变化在时空频域里“放大”给你看。这篇文章我就用Python从零把它写一遍只依赖OpenCV、NumPy和SciPy不套用现成的EVM框架同时把我在调参和复现过程中踩过的坑一并说清楚。这篇内容适合两类人一是刚接触视频处理想搞明白“怎么把一个视频信号拆成时域频域再拼回来”的初学者二是已经在做心率检测、呼吸检测、结构微振动测量之类的工程向开发者想拿Python快速验证一下EVM算法在自己的场景里到底有没有用。我会把从原理、代码到参数调优的完整链路都拆开讲每个步骤都给出选择依据不是单纯甩一个脚本给你。1. 一个“放大看不见的运动”的算法EVM为什么值得自己动手实现1.1 欧拉视角和拉格朗日视角的根本区别传统视频运动分析做得最多的是“找对应点”。比如光流法先在第一帧里找到某个角点再去第二帧里寻找它移动到了哪里这个思路叫拉格朗日视角本质上是跟踪质点本身的位移。EVM用的完全是另一个思路我不关心某个点去了哪我只固定住某一个像素坐标位置观察这个位置上的亮度值随时间怎么变化。比如一面墙壁人眼看毫无动静但每个像素的灰度值其实一直在很小幅度地波动波动频率对应着墙后空调压缩机的振动。这种观察方式叫欧拉视角。这就像你在河边看水流拉格朗日视角是追着一片叶子往下游看它漂到哪欧拉视角是站在桥上盯着固定的一块水面看浪花从那里经过时的起伏规律。EVM之所以叫“欧拉”就是因为它处理的不是运动轨迹而是每个固定像素上的时间序列信号。1.2 泰勒展开为什么放大亮度波动就能放大运动EVM最核心的数学基础其实是大一就能看懂的一阶泰勒展开。假设一维图像信号是 f(x)物体发生了一个微小位移 δ(t)实际采集到的图像是 f(x δ(t))。这个式子可以展开成f(x δ(t)) ≈ f(x) δ(t) · ∂f/∂x其中 ∂f/∂x 是图像在该位置的梯度。注意右边第二项由“位移”乘以“梯度”组成。如果我把这一项放大 α 倍再加回到原信号上得到f(x) (1 α) · δ(t) · ∂f/∂x ≈ f(x (1 α) · δ(t))也就是说刚才那个微小的位移 δ(t)最终在图像上被“撑大”成了 (1α)·δ(t)。这就是欧拉视频放大最核心、也最漂亮的地方不需要知道位移具体是多少只需要把亮度变化的时域信号提取出来放大效果等价于放大了运动幅度。不过直接对像素做这种放大是危险的因为图像里不是所有位置都有清晰的梯度。平坦区域几乎没有梯度放大后没有效果而边缘区域虽然有梯度直接放大又容易导致振铃和噪点。所以EVM还要引入第二个核心工具空间金字塔分解。1.3 真正值得动手实现的应用场景EVM不是实验室里的花架子。最经典的Demo就是用摄像头对着人的脸捕捉皮肤颜色随时间变化的规律从而估计心率。人的皮肤血流随心跳变化时肤色会有极其微弱的改变肉眼完全看不出来但EVM可以把这种颜色变化放大让你直接在视频里看到面部肤色周期性变红变暗。同样可以用在呼吸检测上胸腹部起伏幅度通常很小尤其是在长距离监控场景下普通视频里几乎看不见用EVM放大后腹部起落一目了然。工业上还有用它做桥梁、缆绳、机械结构的微振动分析或者从远处通过窗户玻璃的振动判断室内是否有人说话。可以说任何“视觉上静止但物理上在振动”的场景都适合先拿EVM试一遍。自己动手实现EVM还有一个额外好处你会顺带掌握金字塔分解、时域滤波、视频重建这一整套信号处理链路。这套能力放到很多别的项目里——比如图像融合、风格迁移、视频去隔行——都能复用。2. 环境准备和整体的视频放大管线设计2.1 用到的Python依赖与安装注意事项这个项目用到的库不多就三个opencv-python负责视频读写、图像缩放、金字塔构建。numpy负责多维数组运算和FFT。scipy负责带通滤波器设计主要用Butterworth滤波。安装直接用pip就行。个人建议在虚拟环境里装避免把系统Python搞乱。pip install opencv-python numpy scipy有一个细节要提醒OpenCV的pyrDown/pyrUp函数不同版本对图像尺寸的舍入规则略有差异。如果你用的是4.x版本读取视频帧之后最好确认一下宽高是否为偶数否则在构建金字塔时可能会遇到尺寸不匹配的报错。最简单的处理方式读取层数时强制保持一致后面代码里我会用到dstsize参数显式指定上采样尺寸避免这类问题。2.2 整体管线从视频到放大后的视频EVM的完整流程可以拆成五个阶段读取视频帧把彩色图转成灰度图并归一化到[0,1]区间的浮点数。对每一帧构建拉普拉斯金字塔得到一个多分辨率的图像序列。把视频看作时间和空间两个维度对金字塔每一层沿着时间轴做带通滤波。把滤波后的信号乘以放大系数α加回到原始金字塔系数上获得增强后的金字塔序列。用增强金字塔逐帧重建图像写出视频。这里最关键的是第三和第四步。拉普拉斯金字塔把图像拆成多个空间频带每个频带在时间方向上的信号又经过了带通滤波就相当于同时做了一个时空联合的频域筛选。只有落在你感兴趣的频率范围内的运动或颜色变化才会被放大不会把全频段的噪声全部放出来。2.3 为什么用拉普拉斯金字塔而不是光流法我最早也想用光流法做运动放大但试下来发现有三个问题很难解决。光流法需要计算像素位移位移非常微小的时候估计误差非常大甚至不如直接用泰勒近似来得稳。光流法是逐帧计算的计算量远大于金字塔分解而且对图像噪声极其敏感。光流法得到的是稠密位移场还需要额外做平滑处理流程复杂而金字塔分解天然具备多尺度平滑能力。拉普拉斯金字塔则不同。它的每一层对应着不同的空间频率段居高层捕捉细节和高频梯度居底层捕捉低频平滑结构。微小运动的能量主要分布在梯度明显的边缘处而金字塔分解恰好可以把这些边缘按尺度拆开分别做时域滤波。这样放大后的效果更均匀也不容易出现光流法常见的块状伪影。简而言之拉普拉斯金字塔是EVM算法在空间方向上最重要的支撑结构。没有它直接放大像素时间序列信噪比会差得让人崩溃。3. 核心代码实现逐段拆解3.1 视频读取与帧归一化写代码的第一步是读视频。这里我使用OpenCV的VideoCapture逐帧读取并且在读取时就转成灰度浮点图。为什么要转灰度因为本文主要演示运动放大的核心链路灰度已经足够颜色放大也可以做但需要对每一个颜色通道分别做同样操作后面会提到。import cv2 import numpy as np def read_gray_frames(video_path): cap cv2.VideoCapture(video_path) frames [] while True: ret, frame cap.read() if not ret: break gray cv2.cvtColor(frame, cv2.COLOR_BGR2GRAY) gray gray.astype(np.float32) / 255.0 frames.append(gray) fps cap.get(cv2.CAP_PROP_FPS) cap.release() if not frames: raise ValueError(没能从视频里读取到任何帧检查路径或视频编码。) return frames, fps这里有个容易被新手忽略的点视频帧一般是无符号8位整数范围0~255如果直接做浮点运算很容易产生截断误差和溢出。所以我先转成float32再除以255归一化到0~1否则后面金字塔重建后像素值可能超出8位范围。3.2 构建拉普拉斯金字塔的实现细节拉普拉斯金字塔的构建流程是先构建高斯金字塔再做差分。高斯金字塔每一层都是对上一层做高斯模糊和下采样拉普拉斯金字塔则保存“当前层高斯图”与“上一层上采样到当前尺寸后的图”之间的差异。def build_gaussian_pyramid(img, levels): gp [img] cur img for _ in range(levels - 1): cur cv2.pyrDown(cur) gp.append(cur) return gp def build_laplacian_pyramid(img, levels): gp build_gaussian_pyramid(img, levels) lp [] for i in range(levels - 1): prev_shape (gp[i].shape[1], gp[i].shape[0]) up cv2.pyrUp(gp[i 1], dstsizeprev_shape) lp.append(gp[i] - up) lp.append(gp[-1]) return lp注意dstsize参数必须写成(width, height)不能写成(height, width)。这是OpenCV里比较容易踩坑的一个点。另外最顶层我们放的是高斯金字塔的最后一层也就是一个分辨率极低的低频近似图这样可以保证重建时信息不丢失。这种完整的拉普拉斯金字塔能完美保留原图信息reconstruct_from_laplacian_pyramid是它的逆操作def reconstruct_from_laplacian_pyramid(lp): cur lp[-1] for i in range(len(lp) - 2, -1, -1): dst_size (lp[i].shape[1], lp[i].shape[0]) cur cv2.pyrUp(cur, dstsizedst_size) lp[i] return cur你可能已经发现拉普拉斯金字塔本质上就是一种无损的图像分解方法把图像拆成一个低频基底和若干层差分细节重建时原封不动还原。EVM算法之所以选择它是因为它能让我们在“差分细节层”和“低频基底层”上分别做时域处理灵活度极高。3.3 时域带通滤波如何对金字塔系数做滤波EVM的滤波对象不是整张图像而是金字塔每一层的像素时间序列。也就是说如果视频有N帧金字塔第i层的尺寸是h×w我需要处理N个h×w时间序列每个序列长度是N。最简单的办法是用FFT做理想带通滤波把时间序列变换到频域只保留感兴趣的频带然后逆变换回时域。这样做很好理解但会带来吉布斯效应在频带边界处产生振铃。更好的方案是用SciPy的Butterworth带通滤波器通带更平滑。from scipy.signal import butter, filtfilt def butter_bandpass_filter(data, fps, low, high, order4): nyquist 0.5 * fps low_n low / nyquist high_n high / nyquist b, a butter(order, [low_n, high_n], btypeband) # filtfilt 对整段信号做零相位滤波能有效避免相位偏移 return filtfilt(b, a, data, axis0)这里有几个细节需要说明截止频率必须归一化到奈奎斯特频率也就是帧率的一半。filtfilt是零相位滤波不会让滤波后的信号产生时间偏移。如果换成普通的lfilter波形会整体延时对做生理信号分析影响很大。order一般取4就够了阶数越高滤波器越接近理想带通但数值稳定性越差。如果你的视频帧数太少低于滤波器所需的最低样本数filtfilt会报错。这种情况下建议缩短时间窗或者提高视频帧率。在实际处理时每个金字塔层的数据可以堆叠成一个三维数组维度是[帧数, 高, 宽]再沿第0轴做滤波。这样写起来很直观也方便批量处理。3.4 运动放大与图像重建完整的EVM处理流程我用下面这段代码实现。它把金字塔分解、带通滤波、放大、重建全部串起来。def evm_process(video_path, output_path, levels5, alpha30, low0.5, high3.0): frames, fps read_gray_frames(video_path) n_frames len(frames) h, w frames[0].shape print(f读取 {n_frames} 帧分辨率 {w}x{h}帧率 {fps}) # 第一步逐帧构建拉普拉斯金字塔并按层存储 layer_list [[] for _ in range(levels)] for frame in frames: lp build_laplacian_pyramid(frame, levels) for i in range(levels): layer_list[i].append(lp[i]) layer_data [np.stack(layer_list[i]) for i in range(levels)] # 第二步对每一层做带通滤波并放大 enhanced_layer_data [] for i in range(levels): filtered butter_bandpass_filter(layer_data[i], fps, low, high) enhanced layer_data[i] alpha * filtered enhanced_layer_data.append(enhanced) # 第三步逐帧重建并写出视频 frame_h, frame_w frames[0].shape writer cv2.VideoWriter( output_path, cv2.VideoWriter_fourcc(*mp4v), fps, (frame_w, frame_h) ) for t in range(n_frames): lp_t [enhanced_layer_data[i][t] for i in range(levels)] enhanced_frame reconstruct_from_laplacian_pyramid(lp_t) enhanced_frame np.clip(enhanced_frame, 0.0, 1.0) out_img (enhanced_frame * 255).astype(np.uint8) writer.write(cv2.cvtColor(out_img, cv2.COLOR_GRAY2BGR)) writer.release() print(处理完成输出文件为, output_path)值得注意我在重建前做了np.clip把像素值限制在0到1之间。这个操作非常重要因为放大滤波后叠加后的金字塔系数很可能超出原图范围。如果不做clip最终视频里会出现大量过亮的闪烁点。alpha表示放大倍数。理想情况下alpha越大微小运动放大越明显但噪声也会被一起放大。这个值没有万能设定具体多少合适要看你的视频质量、金字塔层数和拍摄对象我后面会专门讲参数调节策略。现在你可以把任意一段固定机位拍摄的微小运动视频传入evm_process函数试试。如果你手头没有合适素材我建议自己拍一个手机固定在三脚架上对着一个放在桌上、旁边放个蓝牙音箱播放低频正弦波的纸杯录30秒左右你会看到纸杯边缘的微小振动被放大到肉眼清晰可见。4. 参数怎么调才出效果alpha、金字塔层数与频带范围4.1 放大倍数alpha不是越大越好很多第一次跑通EVM的人看到效果不明显第一反应就是把alpha调大。我也这么干过把alpha设到100结果输出视频满屏雪花完全没法看。原因是这样的带通滤波后的信号里不仅有真实微动信号还有传感器噪声和压缩噪声。alpha把真实信号放大的同时也把所有噪声放大。当alpha大到一定程度噪声幅度会盖过真实信号画面就会显得“脏”。一个比较实用的经验对于清晰度较高、噪声较小的视频运动放大的alpha可以取20~50对于压缩严重的视频或手机夜景视频alpha取5~15比较稳妥。如果你只做颜色放大比如心率检测alpha可以用到100以上因为颜色信号的噪声特性与运动信号不同。另外在经典EVM论文的理想化假设里alpha必须满足一个约束放大后的位移不能超过图像结构本身的空间尺度否则泰勒展开近似失效。举个例子一个边缘的梯度变化周期是10个像素你把位移放大了50倍梯度方向上的位移远大于这个周期重建时就会出现严重的鬼影。所以alpha过大的信号往往会看到物体的边缘出现双重轮廓这其实是一个“物理”上的失真。4.2 带通频带的选择先确定你关心的物理现象频率带通频率范围的设置直接决定了算法会放大什么、滤掉什么。这里需要你对自己要观察的物理现象有一个大致频率范围。应用场景频率范围建议说明心率估计0.8 ~ 2.0 Hz对应48~120次/分钟心跳大部分成人心率都落在这个区间呼吸检测0.2 ~ 0.5 Hz对应12~30次/分钟呼吸机械/结构振动根据设备确定比如电机振动通常是5~20 Hz但需要帧率至少达到40以上一般微运动演示0.5 ~ 3.0 Hz比较宽容适合先跑通流程如果你的视频帧率是30fps奈奎斯特频率是15Hz带通上限设到3Hz是完全没问题的。但要注意FFT或Butterworth滤波器的频率分辨率与视频时长有关。30fps下录2秒视频只有60帧能够分辨的频率间隔大约0.5Hz此时设置0.2Hz的带通下界就没有意义因为矩形窗的频率分辨率根本不够。所以录制素材时至少要保证目标频率有5个完整周期以上否则滤波结果不可靠。我在实测中会用这样一个简单办法先把某个关注区域的像素灰度时间序列提取出来做一次FFT看看频谱峰在哪里然后再设置带通范围。这比盲调参数高效得多。4.3 金字塔层数空间尺度与噪声的关系拉普拉斯金字塔层数决定了你把图像分解成多少个空间频带。层数过少不同空间频率的信号混在一起放大效果粗糙层数过多最顶层的分辨率太低重建误差累积反而会让画面变糊。我个人的建议是对于720p及以上分辨率level取5对于480p左右的分辨率level取4如果视频分辨率只有320plevel取3就够了。核心判断标准是金字塔最顶层的最小边长不要小于8像素否则那一层几乎不包含有效结构信息。还有一个多Scale调参技巧不一定所有层都使用同一个alpha。高频层更早的拉普拉斯层次对应细节更容易包含噪声可以适当减小alpha低频层靠后的层可以增大alpha因为它们主要反映整体结构和大幅度运动项。我经常把alpha设置成一个从高频到低频递增的数组比如alphas [10, 15, 25, 35, 35]这样能在保留高分辨率细节的同时避免噪声被过度放大。4.4 一套可直接上手的参数组合参考如果你不想从零开始调可以先照着我这套配置跑然后再根据视频现象慢慢调整。视频类型levelsalphalow(Hz)high(Hz)order人脸心率/肤色变化5800.82.04胸部呼吸运动4400.20.54物体微振动音箱/马达5300.53.04高噪声监控视频4120.32.02需要注意这里给的只是启动值。真实世界的光照条件、压缩格式、运动幅度差异非常大最终参数必须回到“看频谱”和“看输出视频”两个手段上去调。5. 实测中那些文档里不会写的坑5.1 颜色抖动和伪影从哪里来怎么防我跑EVM时最常见的失败现象是画面“呼吸感”特别强整个图像亮度一起一伏但看不到物体边缘放大的微运动。原因通常是带通滤波器通带太宽放大了全局亮度变化而不是局部边缘变化。解决方法是严格限制频带范围同时增加金字塔层数把空间频率划分得更细。还有一种是边缘出现透明拖影看起来像物体在“融化”。这通常是alpha太大导致泰勒展开失效重建时拉普拉斯金字塔的差分系数变得不自然。解决方法是降低alpha或者只对特定金字塔层放大。当画面中出现明显的网格状纹路时多半是金字塔层数太多重建时上采样插值噪声被放大。这时候减少levels并检查cv2.pyrUp的dstsize是否明确设置过。5.2 内存与处理速度的现实问题我的第一版EVM处理代码大概长这样把所有帧、所有层的金字塔系数全部堆到内存里再一次性滤波。这个思路能跑通但极耗内存。如果输入一个720p、300帧的视频level5金字塔系数大概要占掉几百MB内存再乘上中间变量很容易把环境跑崩。优化思路有几种降低输入分辨率。只要不丢失目标微小运动的语义缩放到一半大小往往就够了。滑动窗口处理。不要一次性读入全部帧而是每次处理一个50帧的窗口窗口间重叠10帧避免滤波边界不连续。把处理单元从“整个金字塔层”改成“每个金字塔层的块”分块滤波进一步降低内存峰值。如果追求速度可以用IIR滤波器代替filtfilt每帧只做一次递推更新这样可以做到近似实时。但必须提醒一句filtfilt是零相位滤波数据质量高IIR滤波会引入相位延迟做生理信号提取时需要额外补偿。如果只是做视觉演示IIR完全够用。5.3 怎么验证你放大的是真实信号而不是噪声幻觉EVM输出的视频“看起来有东西在动”不一定说明算法成功了。我见过有人拿一段纯静态噪声视频跑EVM调整alpha后也能看到闪烁但那只是噪声被放大并不是真实微动。验证方法非常简单在放大前后分别取同一个固定像素位置的时间序列画出它们的频谱。如果放大后目标频带对应的峰值显著增高而其他频带变化不大那说明你的EVM确实放大了目标信号。反之如果所有频带都被整体抬高那就是噪声被放大了参数设置有问题。这个检查用Python写非常快import numpy as np def check_spectrum(data, fps): win np.hamming(len(data)) spec np.abs(np.fft.rfft(data * win)) freqs np.fft.rfftfreq(len(data), 1.0 / fps) return freqs, spec把滤波前后的信号丢进去对比一眼就能看出差别。我每次调参数前都会先这样自检一遍能省下大量试错时间。6. 从Demo走向实用的一些扩展思路6.1 实时处理与IIR流式滤波如果想把EVM做到实时金字塔分解和重建在OpenCV里本来就不慢真正拖后腿的是filtfilt这种离线滤波方式。实时场景下可以改用IIR带通滤波器每一帧对金字塔逐层系数做一次迭代更新几乎不增加额外内存。IIR滤波的实现方式是用scipy.signal.lfilter在前一帧滤波器状态的基础上继续往后算。需要注意的是IIR滤波有相位延迟这个延迟会随着滤波器阶数和截止频率变化。如果只做视觉反馈问题不大如果做心跳频段的实时数值输出最好先离线标定一下延迟时间。6.2 从“放大视频”到“直接提取信号”EVM最有价值的用途不一定是要输出放大后的视频而是把视频变成一条可用于分析的时间序列。比如心率检测我不需要真的看到皮肤变红我只需要知道肤色信号在哪个频段存在峰值。这时可以完全跳过金字塔重建过程只取金字塔某一层的某个ROI区域的均值时间序列做FFT找峰即可。这种思路的计算量比完整EVM小一到两个数量级而且对参数不敏感更适合工程落地。我自己做的呼吸检测原型就是直接取金字塔第3层某个区域的灰度均值做一个0.2~0.5Hz的带通滤波然后找峰值间隔准确率已经足够应付大多数场景。6.3 进一步替换为相位放大算法经典EVM基于拉普拉斯金字塔的幅度信息本质上依赖梯度对噪声和大运动的鲁棒性有限。如果你以后需要处理更复杂的场景比如镜头有轻微移动或者目标微动幅度比较大可以考虑“基于相位的运动放大”。这个方向用的是复数可控金字塔或者复数小波变换提取图像的局部相位然后在相位域做滤波和放大。因为相位对幅度变化不敏感所以抗噪性明显更好。不过它的计算复杂度比EVM高不少实现也更复杂。我建议先把经典EVM吃透再决定要不要往相位方向深入。在我自己的项目里EVM算法的代码复杂度不算高但每一个环节的参数都有物理意义这也是我特别喜欢它的原因。把它从零实现一遍相当于把“视频不过是一段三维时空信号”这个抽象概念变成了可以亲手操作的工程经验。刚开始调不出效果很正常建议一次只改一个变量用频谱图去验证每一个假设多试几轮你会慢慢找到感觉。
返回列表