Flowtrace:生物流体中相干结构的简单可视化（2017）
原文链接：https://doi.org/10.1242/jeb.162511
前段时间一直在说要更新一下有关流场的一些函数，这部分内容既有python的，也有MATLAB的，所以就在这篇文章的底部附上对应的一些代码的应用吧！这篇文章是16年的，特别有意思，其主要思路就是根据所拍摄的生物视频的多帧进行流场流线的预测和绘制，主要功能是定性地了解生物流场，类似于粒子图像测速技术（PIV）中的可视化，不过它的优势在于没有PIV复杂，但是同样也不具备PIV地定量测量，所以我们开始论文阅读和后面的一些代码学习吧！
在文中除了对英文的翻译外，我将对应的文献进行了注明，另外也会添加一些我自己的看法，我会用斜体来表示。
1 摘要
我们提出了一种简单、直观的算法，用于可视化时变流场，可以在用户干预最少的情况下揭示复杂的流动结构。我们将这种技术应用于各种生物系统，包括无脊椎动物的游泳水流和成群昆虫的集体运动。我们将我们的结果与实验难度更大、数学上更复杂的识别流体模式的技术进行了比较，并建议我们的工具代表了一个重要的“中间地带”,允许实验人员容易地确定一个系统是否表现出有趣的流动模式和连贯的结构，而不需要求助于更密集的技术。除了信息丰富之外，我们的工具生成的可视化效果通常是引人注目的和优雅的，直接从视频中说明连贯的结构，而不需要计算叠加。我们的工具可以作为完整记录的开源代码在www.flowtrace.org大学的MATLAB、Python或ImageJ上获得。
这是一个非常实用的工具，我是最近在实验生物学杂志（JEB）上找到的，这里提到说明了该软件工具的应用对象是：无脊椎动物以及成群昆虫，后面会有图片举例说明，它的应用所要实现的目标是识别所探测物体所产生的流体模式，通俗点讲就是可视化流体运动场，文章页数不超过十页，还请大家耐心阅读哦。
2 简介
流动可视化是研究许多不同领域(包括生态学、1生物力学、2和分子生物学)中出现的生物过程的一项基本技术。3，4观察许多示踪粒子的轨迹——无论是荧光珠、微泡、移动细胞内的小泡还是阻塞细菌——可以揭示生物运动(如分子运输、纤毛流动或主动力产生)在许多长度和时间范围内的协调方式的丰富信息。5-8此外，流动可视化有助于许多标准技术的表征和二次验证。 如微流体和流式细胞仪9，10。
最简单的定性流动可视化技术包括观察被动标量的运动，如染料或烟雾，因为它们是由流体平流的。 这些实验的优点是执行起来相对简单，并且通常可以直接洞察流动的整体结构和混合特性。11但是，在扩散效应占主导地位的长度和时间范围内，或者以不同过程的时间范围内的大分离为特征的流入，这些研究的结果可能难以解释。12，13因此，在许多情况下，基于示踪粒子运动的定量流动表征技术是优选的。14–16
然而，标准的定量流动可视化技术——粒子跟踪和粒子图像测速术(PIV)——很难实施和优化，因此严格的流体动力学可视化技术在许多实验环境中仍然受到禁止。13虽然PIV和相关技术已广泛应用于某些系统并进行优化，如鱼类游泳17或血流力学的研究，8但在不太成熟的环境中，这些技术需要对仪器进行专门的修改，如构建激光光幕或点探针。4，18因此， 依赖于系统的技术通常是必要的，19特别是当粒子运动在数据中仅部分可见时，例如从图像条纹9或由于有限的景深造成的失焦漂移。20这个问题对于非传统介质中的流动更为突出，例如在岩群和兽群的集体运动中。21–23
然而，流体力学中的许多标准概念——例如涡流、射流和湍流——对于整个科学领域的研究人员来说是众所周知的。 即使不需要定量流动表征工具，他们也可能在其数据中识别出这些特征的可能存在。24，25这表明，对于 基于染料的定性技术不可用的系统，额外的流动可视化工具是必要的，但定量示踪粒子研究是不必要的。
在这里，我们介绍了Flowtrace，一种直观的定性可视化技术，可以揭示实验视频中各种流动结构的存在，使其既可以作为存在/不存在研究的主要分析工具，也可以促进更复杂的流动量化技术的使用。我们的技术完全基于输入数据的图像处理，而不是像涡度那样的标量场的数值重建，这使得它可以用作生物系统的直接“首次通过”表征技术，其中传统技术要么不需要支持定性观察，要么由于所涉及的长度和时间尺度而非常困难。虽然我们的技术简单直观，但我们无法在文献中找到以前的报告——重要的是，我们发现它在研究广泛的生物数据集中的模式方面具有令人惊讶的效用。
更为详细地阐述了Flowtrace工具创建的背景因素是：染料技术的不足以及PIV技术的复杂，意义是生物的模式观察需要可支持的简便工具。
3 方法
3.1 算法
我们的算法概括了一种从视频生成长曝光照片的常见技术，在该技术中，获取时间序列的最大(或最小)强度投影，以便为在暗背景下移动的明亮对象生成路径线，从而在图像上产生运动条纹。5该技术以前曾被用于创建星迹，这是一种流行的天文可视化技术，通过获取稳定的夜空视频的最大强度投影而生成。26
我们对这一技术的主要扩展是拍摄一段长视频，然后对连续的小组帧进行最大强度投影，以便生成显示不同时间路径线的序列图像。得到的路径线图像的时间序列显示了路径线的形状如何随时间变化。具体地说，对于100帧视频，具有30帧轨迹的流迹视频由帧1-31、2-32、3-33等的最大强度投影组成。当投影第70帧到第100帧时，序列结束，从而生成70帧视频。除了此基本运算外，还可以将各种其他运算与最大强度投影运算相结合，以产生更好的结果。图1C以图形方式说明了该过程。
象征性地，让电影的每一帧是像素值和位置的向量vij[t]，其中t = 1，2，…，N表示由N帧组成的视频中一帧的索引，i和j表示图像中一个像素的坐标。假设示踪粒子是在黑暗背景下运动的亮色物体。在这种情况下，我们根据前向卷积算子定义一系列最大强度投影p[t],
其中M < N是视频中帧的某个子集。随着时间索引t向前“滑过”连续的索引1、2、3，…，最大强度投影是在M帧的连续运行中获得的，每个帧相差两个图像(第一个和最后一个)。这导致一组最大强度投影，ptij，t = 1，2，…，N-M，它们构成了从原始数据集生成的新视频。重要的是，生成的视频中的帧数(N M)几乎等于原始视频中的帧数(N)。在本文中，我们将M幅图像的每个子序列称为子堆栈，将单个粒子在M帧间移动所占据的位置序列称为路径线。参数M是粒子轨迹线在每一帧中可视化的时间刻度，它是用户使用该工具必须指定的唯一参数。
在卷积期间，投影算子P“滑过”整个帧序列，在重叠的M帧组中对视频进行操作。该操作符可以与其他预处理操作组合，以达到不同的效果；在下面描述的代码中，定义的其他操作包括中值减法(移除缓慢移动的对象)、颜色反转(针对在较亮背景下移动的较暗对象)、差分加权(对M帧序列中的每一帧进行不同量的着色或变暗，以便显示指示时间的路径线上的梯度)，以及成对差分(隔离移动速度超过1 px/帧的对象)。这些操作中的每一个都象征性地表示卷积之前的合成，使得最终的图像序列是((P◦G)*v)[t]，其中G是预处理操作。
四个大法宝：中值减法；颜色反转；差分加权；成对差分。
3.2 软件包和选项
我们的算法是作为“Flowtrace”实现的，这是www.flowtrace.org大学的一个用于MATLAB、Python或ImageJ的开源包。那里提供了完整的教程和示例图像集。表1总结了可用于代码的主要用户指定的参数和选项；可选参数在MATLAB中作为struct对象传递，在Python中作为关键字参数传递，在ImageJ的GUI中作为复选框传递。（表1看不清楚可以下载原文查看啊）
通常一个数据集有两个截然不同的速度标度 ；一个用于示踪粒子，另一个用于逐渐漂移、大流量等的慢分量。在这种情况下，当投影长度M设置为使快速粒子在视场内移动很远，而以较慢速度移动的对象移动相对较少时，Flowtrace的性能最佳。这对于图1中所示的两种产生喂食电流的生物体的轨迹线来说是正确的，这些轨迹线是针对投影周期绘制的，该投影周期与示踪粒子的通过时间相比是长的，但是与每个生物体身体的渐进运动相比是短的。结果，每个生物体的解剖特征在图像中保持清晰。
在某些情况下，这些时间尺度没有很好地分开 (导致慢速移动物体的运动模糊)，或者图像中存在静止物体和障碍物，使得路径模糊。对于这些情况，在执行投影之前对每个子堆栈应用背景减法操作是有用的。对于相对于示踪粒子缓慢移动但足够快以显示明显的运动模糊的对象，最积极的背景选项是“获取差异”，它在应用投影之前获取所有连续图像之间的成对差异。然而，如果一些示踪粒子相对于其他粒子移动缓慢，这些粒子也会从图像中消失。对于几乎静止的背景对象，“减法优先”或“减法中值”基于对象的类型产生类似的结果。有时，在生成的路径中突出显示时间的方向性是很方便的，特别是在使用输出时间序列的静止帧进行分析时。在这种情况下，可以通过使用可选参数“颜色系列”应用跨时间的颜色渐变，或者使用可选参数“渐变尾部”应用线性强度渐变来指示方向。
Flowtrace的主要涉及的四种方法：
- 背景减法与获取差异：解决运动缓慢物体及背景障碍物所产生的图像模糊问题；
- 颜色渐变：将时间序列的变化显示在静止的帧数上；
- 减法中值：凸显多种物体中低速物体的运动轨迹。
4. 结果和讨论
我们已经将Flowtrace应用于各种生物系统，并且它被证明在解释各种情况下出现的复杂的不稳定流动动力学方面是惊人地有效。通过改变投影时间间隔τ(以及路径长度),我们已经研究了在大范围的长度和时间尺度上发生的生物现象。下面简要描述的每个实验的方法将在补充材料中进一步讨论。
图1A显示了海星幼虫流迹电影的三个代表性帧，在之前的工作中，我们已经展示了当海星幼虫不断调整其喂食电流时，在其身体周围产生的动态涡流阵列。27生成视频的完整流迹电影显示了随着时间的推移，不同的涡流模式平滑地过渡到另一个，表明这些不同的动态流模式如何由于动物的控制机制而演变(见参考文件的补充视频11)。27在图1B中，我们探索了原生动物斯藤托耳种在过滤器喂食过程中产生的类似纤毛流。(尺寸50米)。视频和图像(τ = 3秒)捕捉到了生物缓慢旋转其茎时藻类颗粒的螺旋运动(补充视频S1)。
在海葵Aiptasia pallida (∼ 1 mm)中，一个反转彩色视频(τ = 4分钟)显示了当动物蠕动地将含有示踪珠的水泵入其体腔时水射流的破裂(图2C，补充视频S3)。在一个类似的月球蜗牛面盘虫(1毫米)的视频中，Flowtrace通过分别投影每个通道来处理DSLR相机生成的彩色视频，生成由游泳动物创建的双极流场形成的真彩色视频(图2D，补充视频S4)。
除了被动流体示踪粒子，Flowtrace还可应用于主动粒子和生态数据。在一个爬行的军蚁群中(∼ 5 mm，用手持iPhone在田间拍摄)，Flowtrace显示了蚂蚁路径的收敛(τ = 1 s，图2B，补充视频S2)。类似地，在一群宽1米的飞蠓中，应用可选参数颜色系列显示沿蠓路径的颜色梯度(代表时间),允许可视化收紧封锁(图2A，补充视频S1，原始数据取自参考论文)。28
Flowtrace使用的简单滑动投影技术在生物科学和流体动力学文献中似乎鲜为人知，尽管它可以很容易地实现。Flowtrace可以再现几项研究的核心定性结论，包括我们最近对幼海星游泳的研究，27，29，其中动物产生的不同进食和游泳涡流阵列的关键观察首先基于Flowtrace视频确定。我们已经发现，我们的方法对于馈送电流的研究特别有效，因为在平流和场变化的时间尺度之间有很大的间隔，使得示踪粒子在流场经历进一步变化之前有足够的时间绘制出流场的结构。此外，喂食现象通常涉及小的长度尺度和长的时间尺度，传统的染料平流可视化技术将由于快速扩散混合而失败。
我们将Flowtrace与其他识别流体结构的技术进行了比较，我们发现Flowtrace可以定性地再现使用涡量等值线或有限时间李亚普诺夫指数的更复杂动态分析的结果30，31(见补充分析)。我们的工具可用于 Fiji/ImageJ、Python 2和3以及MATLAB的高效多线程实现。我们设想一个广泛的社区——从显微镜学家到生态学家到流体物理学家——会发现Flowtrace是有用的，所以完整的源代码和文档可以在http://www.flowtrace.org获得，或者在GitHub上请求。
还有几幅图文中没有着重点到，所以也就不放在这里了！最后把一些有关*.flo**后缀的流场速度文件的可视化代码给介绍一下。*
首先我们需要注意.flo文件中坐标的定义，见下图！和我们数组的定义非常的相近，向右为水平正值，向下为竖直正值。
MATLAB的readFlowFile和flowToColor
参考链接：https://vision.middlebury.edu/flow/
readFlowFile的源代码 ： 读取.flo文件为三维的可处理数组
function img = readFlowFile(filename)
% readFlowFile read a flow file FILENAME into 2-band image IMG
if isempty(filename) == 1
error('readFlowFile: empty filename');
end;
idx = findstr(filename, '.');
idx = idx(end);
if length(filename(idx:end)) == 1
error('readFlowFile: extension required in filename %s', filename);
end;
if strcmp(filename(idx:end), '.flo') ~= 1
error('readFlowFile: filename %s should have extension ''.flo''', filename);
end;
fid = fopen(filename, 'r');
if (fid < 0)
error('readFlowFile: could not open %s', filename);
end;
tag = fread(fid, 1, 'float32');
width = fread(fid, 1, 'int32');
height = fread(fid, 1, 'int32');
% sanity check
if (tag ~= TAG_FLOAT)
error('readFlowFile(%s): wrong tag (possibly due to big-endian machine?)', filename);
end;
if (width < 1 || width > 99999)
error('readFlowFile(%s): illegal width %d', filename, width);
end;
if (height < 1 || height > 99999)
error('readFlowFile(%s): illegal height %d', filename, height);
end;
nBands = 2;
% arrange into matrix form
tmp = fread(fid, inf, 'float32');
tmp = reshape(tmp, [width*nBands, height]);
tmp = tmp';
img(:,:,1) = tmp(:, (1:width)*nBands-1);
img(:,:,2) = tmp(:, (1:width)*nBands);
fclose(fid);
flowToColoe的源代码 ： 将readFlowFile读取后的数组变为伪色图，该转换标准可以自行查阅17年的这篇论文（包含较多光流方面的知识非常不错哦！）： Scene Flow Estimation : A Survey.
function img = flowToColor(flow, varargin)
% flowToColor(flow, maxFlow) flowToColor color codes flow field, normalize
% based on specified value,
%
% flowToColor(flow) flowToColor color codes flow field, normalize
% based on maximum flow present otherwise
UNKNOWN_FLOW_THRESH = 1e9;
UNKNOWN_FLOW = 1e10; %
[height widht nBands] = size(flow);
if nBands ~= 2
error('flowToColor: image must have two bands');
end;
u = flow(:,:,1);
v = flow(:,:,2);
maxu = -999;
maxv = -999;
minu = 999;
minv = 999;
maxrad = -1;
% fix unknown flow
idxUnknown = (abs(u)> UNKNOWN_FLOW_THRESH) | (abs(v)> UNKNOWN_FLOW_THRESH) ;
u(idxUnknown) = 0;
v(idxUnknown) = 0;
maxu = max(maxu, max(u(:)));
minu = min(minu, min(u(:)));
maxv = max(maxv, max(v(:)));
minv = min(minv, min(v(:)));
rad = sqrt(u.^2+v.^2);
maxrad = max(maxrad, max(rad(:)));
fprintf('max flow: %.4f flow range: u = %.3f .. %.3f; v = %.3f .. %.3f\n', maxrad, minu, maxu, minv, maxv);
if isempty(varargin) ==0
maxFlow = varargin{1};
if maxFlow > 0
maxrad = maxFlow;
end;
end;
u = u/(maxrad+eps);
v = v/(maxrad+eps);
% compute color
img = computeColor(u, v);
% unknown flow
IDX = repmat(idxUnknown, [1 1 3]);
img(IDX) = 0;
综合示例：
addpath functionsLib;%存放上面两个函数的文件夹
clc
clear
flow_1 = readFlowFile('文件名.flo');
img_1 = flowToColor(flow_1);
imshow(img_1);%图4所示
python的cv2和numpy
当然python版本的上面两个函数也是有的不过我们来点不一样的，cv2这个库里面有readOpticalFlow文件很好用，至于颜色转换的函数可以使用numpy的平方根sqrt计算，最后使用matplotlib.pyplot画图图像。输出结果就不展示了，喜欢的朋友自己搞一搞。
看代码：
import cv2
import numpy as np
import matplotlib.pyplot as plt
Flows = cv2.readOpticalFlow(flows)#flows为.flo文件的路径
Flow = np.sqrt(Flow[:,:,0]**2+Flow_z[:,:,1]**2)#注意用sqrt哦！
plt.matshow(Flow, cmap=plt.cm.viridis)#viridis看的最舒服了
plt.show()# 这个不要忘了
25
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
