1. 项目概述从时序信号到图像识别的故障诊断新思路格拉姆角场Gramian Angular Field GAF结合轴承故障诊断这个组合在工业预测性维护领域已经不算新鲜但对于刚接触的同学来说看到东南大学相关的代码和数据集第一反应可能还是有点懵好好的振动信号为什么要费劲转换成图像直接用深度学习模型处理一维时序数据不行吗我最初也有这个疑问。直到在实际项目中面对来自不同工况、带有强噪声的轴承振动数据传统时频分析方法比如FFT、小波变换的特征提取稳定性遇到了瓶颈而基于GAF的方法展现出了独特的优势。简单来说GAF的核心思想是将一维时间序列通过坐标变换映射到极坐标系再通过三角运算构造出一种类图像矩阵表示。这种表示方法巧妙地将时间序列的时序依赖性和数值关系“凝固”在一张图上使得后续可以借助在图像识别领域非常成熟的卷积神经网络CNN来进行特征学习和分类。这对于轴承故障诊断而言相当于开辟了一条“降维打击”的新路径——我们不再需要手工设计复杂的时域、频域、时频域特征而是让CNN从这种“翻译”过来的图像中自动学习故障的视觉模式。东南大学在机械故障诊断领域的研究一直走在前列其公开的代码和数据集为初学者和研究者提供了极佳的学习范本。本次解读聚焦于“数据集解读”部分因为这是整个流程的基石。如果数据都理解错了后面的模型构建、训练调参都是空中楼阁。我们将深入拆解代码中数据加载、预处理、以及最关键的一步——如何将原始的振动信号样本转换为GAF图像的全过程并分享我在复现和扩展过程中踩过的坑和总结的经验。2. 核心思路与方案选型为什么是GAF在深入代码之前我们必须先搞清楚方案选型背后的逻辑。轴承故障诊断本质上是一个模式识别问题我们需要从传感器采集的振动信号中区分出“正常”、“内圈故障”、“外圈故障”、“滚动体故障”等不同状态。传统方法流程固定原始信号 - 数字滤波降噪 - 特征提取如均方根、峭度、频谱峰值 - 特征选择/降维 - 输入分类器如SVM、随机森林。这个流程的瓶颈在于“特征提取”环节。手工设计的特征严重依赖专家经验且对于变工况、变负载、强噪声的场景泛化能力往往不足。深度学习提供了一种端到端的解决方案但直接将一维振动信号喂给1D-CNN或RNN有时难以充分捕捉复杂的时序动态和长期依赖。这时GAF的优势就体现出来了。它的转换过程可以概括为两个核心步骤归一化与极坐标映射将一维时间序列的数值归一化到[-1, 1]或[0, 1]区间然后将每个数据点视为在单位圆上的一个点其角度由归一化后的值通过反余弦函数决定半径固定为1或由时间戳决定。这一步将时序信息编码到了极角中。生成格拉姆矩阵通过计算每两个点之间的三角和或差的余弦值生成一个格拉姆矩阵。这个矩阵是一个对称矩阵其元素反映了原始序列中任意两点之间的时序关系。这个矩阵就可以被视作一张灰度图像。为什么选择这个方案保留时序信息与直接将序列排列成图像不同GAF的转换过程本质上是时序相关的矩阵中的每个点都包含了两个原始时间点之间的关系。适合CNN处理生成的GAF图像是结构化的、局域相关的二维数据这与CNN擅长的处理对象如图像完美契合。CNN可以高效地从中提取空间层次化特征。对幅度缩放具有不变性由于先进行了归一化GAF对信号的整体幅度变化不敏感更关注信号的形状和相对变化这在工业环境中非常实用因为设备负载变化会导致信号幅度整体漂移。在东南大学的代码实现中通常采用GAF的两种变体格拉姆角和场GASF和格拉姆角差场GADF。简单理解GASF使用余弦和图像更强调序列的整体趋势GADF使用正弦差对序列的局部变化和梯度更敏感。代码中往往会同时生成这两种图像或者选择其中一种作为输入这需要根据具体数据特性进行实验。3. 数据集深度解读与预处理实战拿到一个故障诊断数据集绝不能直接扔进模型。正确的打开方式是先像侦探一样审视它。东南大学常用的数据集包括经典的CWRU凯斯西储大学轴承数据中心的数据也可能是其自有实验台的数据。我们以CWRU数据集为例进行深度解读因为它的结构清晰应用广泛。3.1 数据集结构与物理意义剖析CWRU数据集的目录结构通常按驱动端风扇端、故障直径、负载工况来组织。例如CWRU/ ├── 12k Drive End Bearing Fault Data/ # 12kHz采样驱动端数据 │ ├── Ball007/ # 滚动体故障直径0.007英寸 │ ├── IR007/ # 内圈故障 │ ├── OR007/ # 外圈故障 │ └── Normal/ # 正常状态 └── 48k Drive End Bearing Fault Data/ # 48kHz采样数据每个子文件夹里是多个.mat文件每个文件对应一次采样记录通常包含一个名为DE驱动端加速度的变量也可能包含FE风扇端和BA基座加速度数据。关键参数解读采样频率Fs常见12kHz和48kHz。这决定了信号的最高分析频率根据奈奎斯特定理为Fs/2。对于轴承故障特征频率通常几百Hz到几千Hz12kHz通常足够。故障直径如0.007、0.014、0.021英寸。故障尺寸直接影响振动信号的冲击强度和调制现象。负载如0HP、1HP、2HP、3HP。电机负载不同轴承的受力状态不同故障特征频率的幅值会受负载调制这是模型泛化能力的重要考验。代码中的数据加载环节通常使用scipy.io.loadmat来读取.mat文件。这里第一个注意事项就来了一定要确认加载后数据的维度和变量名。有时数据会被多层嵌套需要用.item()或索引来取出真正的振动信号数组。import numpy as np from scipy.io import loadmat # 示例加载一个.mat文件 file_path ‘path/to/your/data/Normal_0.mat‘ mat_data loadmat(file_path) # 关键查看mat文件中所有变量名 print(mat_data.keys()) # 通常振动数据存储在 ‘DE‘ 这个键下 vibration_signal mat_data[‘DE‘].flatten() # 确保是一维数组 print(f“信号长度{len(vibration_signal)} 采样频率假设为 12kHz“)3.2 数据切片与样本构建策略原始数据文件往往很长比如12kHz采样下10秒就是12万个点我们需要将其切割成多个固定长度的样本用于训练和测试。样本长度segment_length的选择是一个需要权衡的参数太短如1024点可能无法包含一个完整的故障冲击周期信息不充分。太长如8192点或更长计算GAF图像时矩阵尺寸过大n x n计算和存储开销剧增且可能包含过多的冗余信息或多种状态的混合。经验值对于CWRU的12kHz数据我通常选择2048或4096个点作为一个样本。这对应约0.17秒或0.34秒的数据足以捕捉到几次故障冲击同时矩阵尺寸可控2048x2048的图像已经需要约33MB内存存储为float64通常需要下采样或使用更小的切片。代码中的切片操作需要注意重叠问题。为了增加样本数量可以采用重叠切片。例如步长stride设为segment_length//2即50%的重叠率。def create_samples(signal, segment_length2048, stride1024): “”“将长序列切割成固定长度的样本。 Args: signal: 一维振动信号数组。 segment_length: 每个样本的长度。 stride: 滑动步长。 Returns: samples: 形状为 (n_samples, segment_length) 的二维数组。 ”“” n_samples (len(signal) - segment_length) // stride 1 samples np.zeros((n_samples, segment_length)) for i in range(n_samples): start i * stride end start segment_length samples[i] signal[start:end] return samples注意务必确保每个样本的标签是正确的。如果从一个“内圈故障”的数据文件中切出100个样本那么这100个样本的标签都应该是“内圈故障”。在构建最终数据集时需要将不同故障类型、不同工况的样本和标签分别堆叠起来并记得打乱顺序在划分训练集和测试集之后打乱而不是之前。3.3 数据标准化被忽视的关键一步在将样本送入GAF转换之前对每个样本进行独立的标准化至关重要。这是因为GAF的第一步——归一化到[-1,1]区间对输入数据的尺度非常敏感。如果不同样本的绝对幅值差异很大这在变工况数据中很常见直接使用全局的归一化参数会导致部分样本的信息被压缩。正确的做法是对每个样本进行局部标准化常用方法是减去均值除以标准差Z-score标准化或者最小-最大归一化到[-1,1]。这样能保证每个样本自身都被规范到相同的尺度突出了其内部的相对变化模式这正是GAF想要捕捉的。def normalize_sample(sample, method‘zscore‘): “”“标准化单个样本。 Args: sample: 一维数组一个振动信号样本。 method: ‘zscore‘ 或 ‘minmax‘。 ”“” if method ‘zscore‘: mean np.mean(sample) std np.std(sample) if std 1e-10: # 防止除零 std 1.0 return (sample - mean) / std elif method ‘minmax‘: min_val, max_val np.min(sample), np.max(sample) if max_val - min_val 1e-10: return sample * 0 return 2 * (sample - min_val) / (max_val - min_val) - 1 else: raise ValueError(“Method must be ‘zscore‘ or ‘minmax‘“)4. GAF图像生成核心代码逐行解读这是整个流程的技术核心。我们将结合代码详细解释每一步的数学含义和实现细节。4.1 极坐标映射从数值到角度假设我们有一个已经标准化到[-1, 1]区间的样本X {x1, x2, ..., xn}。GAF映射的第一步是计算每个点对应的角度φ。公式为φ_i arccos(x_i), 其中x_i ∈ [-1, 1] 因此φ_i ∈ [0, π]。为什么用反余弦因为它是一个在[-1,1]区间上单调递减的函数能将数值唯一地映射到[0, π]的角度空间。同时由于余弦函数在[0, π]上是单调的这个映射是可逆的。import numpy as np def to_polar_coordinates(normalized_sample): “”“将归一化后的样本转换为极坐标角度。 Args: normalized_sample: 归一化到[-1,1]的一维数组。 Returns: phi: 对应的角度数组范围[0, pi]。 ”“” # 防止数值误差导致归一化值略微超出[-1,1]范围 normalized_sample np.clip(normalized_sample, -1, 1) phi np.arccos(normalized_sample) return phi4.2 生成格拉姆矩阵GASF与GADF得到角度数组φ后我们计算格拉姆矩阵。这里以**格拉姆角和场GASF**为例其元素定义为GASF_ij cos(φ_i φ_j)这个定义可以展开为cos(φ_i)cos(φ_j) - sin(φ_i)sin(φ_j)。注意cos(φ_i)就是我们的原始归一化值x_i。因此GASF矩阵可以直接通过原始归一化数据计算无需显式计算角度效率更高GASF X^T · X - sqrt(1 - X^2)^T · sqrt(1 - X^2)其中X是归一化样本向量sqrt(1 - X^2)即sin(φ)。格拉姆角差场GADF的定义为GADF_ij sin(φ_i - φ_j)其高效计算方式为GADF sqrt(1 - X^2)^T · X - X^T · sqrt(1 - X^2)def gramian_angular_field(sample, method‘sum‘, scaleNone): “”“计算样本的格拉姆角场。 Args: sample: 一维数组一个振动信号样本建议已归一化。 method: ‘sum‘ 对应 GASF ‘difference‘ 对应 GADF。 scale: 是否缩放图像到[0, 255]区间用于可视化。默认不缩放。 Returns: GAF: 二维矩阵即GAF图像。 ”“” # 确保输入是一维且为浮点型 sample sample.astype(np.float64).flatten() n len(sample) # 方法1通过角度计算直观但稍慢 # phi np.arccos(np.clip(sample, -1, 1)) # if method ‘sum‘: # # GASF # cos_sum np.cos(np.add.outer(phi, phi)) # return cos_sum # else: # # GADF # sin_diff np.sin(np.subtract.outer(phi, phi)) # return sin_diff # 方法2通过三角恒等式高效计算推荐 sample_clipped np.clip(sample, -1, 1) # 计算 sin(phi) sin_phi np.sqrt(1 - sample_clipped ** 2) # 注意这里隐含假设phi在[0, pi]sin(phi)0 if method ‘sum‘: # GASF cos(phi_i phi_j) x_i * x_j - sqrt(1-x_i^2)*sqrt(1-x_j^2) outer_prod np.outer(sample_clipped, sample_clipped) sin_outer np.outer(sin_phi, sin_phi) gaf outer_prod - sin_outer elif method ‘difference‘: # GADF sin(phi_i - phi_j) sqrt(1-x_j^2)*x_i - x_j*sqrt(1-x_i^2) # 利用外积计算 sin_phi^T * X - X^T * sin_phi term1 np.outer(sin_phi, sample_clipped) term2 np.outer(sample_clipped, sin_phi) gaf term1 - term2 else: raise ValueError(“Method must be ‘sum‘ or ‘difference‘“) if scale is not None: # 将GAF值线性缩放到[0, scale]区间便于保存为图像 gaf_min, gaf_max gaf.min(), gaf.max() if gaf_max - gaf_min 1e-10: gaf (gaf - gaf_min) * scale / (gaf_max - gaf_min) else: gaf np.zeros_like(gaf) return gaf4.3 图像下采样与存储优化直接生成2048x2048的GAF图像对于大批量训练来说是巨大的内存和计算负担。一个常见的优化技巧是在生成GAF之前先对一维样本进行下采样。例如将2048点的样本通过滑动平均或直接每隔N个点采样的方式降到更短的序列长度M如64 128 256。这样生成的GAF图像尺寸就是M x M大大减少了计算量。from scipy import signal def downsample_sample(sample, target_length256): “”“下采样样本到目标长度。 Args: sample: 一维样本。 target_length: 下采样后的目标长度。 Returns: 下采样后的一维数组。 ”“” original_length len(sample) if original_length target_length: return sample # 方法1简单的线性插值重采样 # return np.interp(np.linspace(0, original_length-1, target_length), # np.arange(original_length), sample) # 方法2使用scipy.signal.resample基于FFT更适用于带限信号 return signal.resample(sample, target_length)存储格式生成的GAF矩阵是浮点型二维数组。如果直接保存为.npy文件体积较大。通常有两种处理方式实时生成在模型训练的数据加载器DataLoader中实时将一批一维样本转换为GAF图像。这种方式灵活不占用额外磁盘空间但会增加每个epoch的训练时间。预处理保存将所有训练集和测试集的样本预先转换为GAF图像并保存为图像文件如.png.jpg或压缩的数组文件如.npz。这种方式训练速度快但需要大量的磁盘空间。对于256x256的灰度图保存为uint8的PNG格式每张图约65KB1万张图就是650MB尚可接受。在东南大学的代码中为了流程清晰和实验可复现通常采用第二种方式即先预处理生成一个完整的图像数据集。5. 完整数据处理管道构建与经验分享将上述所有步骤串联起来就构成了从原始.mat文件到最终GAF图像数据集的数据处理管道Pipeline。这个管道的健壮性和效率直接影响后续实验的顺利进行。5.1 管道构建示例下面是一个简化的端到端管道示例它遍历指定目录下的所有.mat文件生成对应标签的GAF图像并保存到以故障类别命名的文件夹中。import os import numpy as np from scipy.io import loadmat import cv2 # 用于保存图像 import argparse def build_gaf_dataset(data_root_dir, output_img_dir, segment_len2048, stride1024, downsampled_len256, gaf_method‘sum‘): “”“构建GAF图像数据集。 Args: data_root_dir: 原始.mat数据根目录子文件夹为不同故障类别。 output_img_dir: 输出图像目录内部会按类别创建子文件夹。 segment_len: 样本切片长度。 stride: 切片步长。 downsampled_len: 下采样目标长度。 gaf_method: ‘sum‘ (GASF) 或 ‘difference‘ (GADF)。 ”“” # 获取故障类别子文件夹名 fault_classes [d for d in os.listdir(data_root_dir) if os.path.isdir(os.path.join(data_root_dir, d))] fault_classes.sort() print(f“发现故障类别{fault_classes}“) for class_idx, class_name in enumerate(fault_classes): class_dir os.path.join(data_root_dir, class_name) output_class_dir os.path.join(output_img_dir, class_name) os.makedirs(output_class_dir, exist_okTrue) mat_files [f for f in os.listdir(class_dir) if f.endswith(‘.mat‘)] print(f“处理类别 ‘{class_name}‘ 共 {len(mat_files)} 个文件。“) sample_count 0 for mat_file in mat_files: file_path os.path.join(class_dir, mat_file) try: mat_data loadmat(file_path) # **关键确认数据键名这里假设为 ‘DE‘** vibration_signal mat_data[‘DE‘].flatten() except Exception as e: print(f“加载文件 {file_path} 失败{e}“) continue # 1. 创建样本切片 samples create_samples(vibration_signal, segment_len, stride) # 2. 对每个样本进行处理 for i, sample in enumerate(samples): # 2.1 样本标准化 (Z-score) norm_sample normalize_sample(sample, method‘zscore‘) # 2.2 下采样 (可选但强烈推荐) if downsampled_len and len(norm_sample) downsampled_len: norm_sample downsample_sample(norm_sample, target_lengthdownsampled_len) # 2.3 生成GAF图像 gaf_matrix gramian_angular_field(norm_sample, methodgaf_method) # 2.4 缩放到[0, 255]并转换为uint8以便保存为图像 gaf_img ((gaf_matrix 1) / 2 * 255).astype(np.uint8) # 假设GAF值在[-1,1] # 2.5 保存图像文件名包含类别和索引信息 img_filename f“{class_name}_{sample_count:06d}.png“ img_path os.path.join(output_class_dir, img_filename) cv2.imwrite(img_path, gaf_img) sample_count 1 print(f“ 类别 ‘{class_name}‘ 处理完成生成 {sample_count} 张图像。“)5.2 实操心得与避坑指南数据平衡问题不同故障类别的原始数据文件可能数量不等切片后样本数差异可能更大。务必在构建最终数据集时检查各类别的样本数量。如果严重不平衡需要考虑过采样如对少数类样本进行随机滑动窗口的多次切片、欠采样或使用类别权重。GAF图像的可视化检查在管道运行后一定要随机抽查几张生成的GAF图像用matplotlib显示出来看看。正常的、不同故障的GAF图像应该呈现出不同的纹理模式。如果所有图像看起来都是模糊一片或噪声很可能是在数据标准化或GAF计算环节出了问题。下采样参数的权衡downsampled_len是核心超参数。太小如32会丢失过多细节太大如512则计算成本高。建议从128或256开始尝试。下采样方法也会影响效果scipy.signal.resample比简单的线性插值更能保留频域信息。内存管理处理大规模数据时避免一次性将所有数据加载到内存。上述管道是“流式”处理处理一个文件保存一批图像是更安全的方式。如果使用实时生成策略在DataLoader中要确保转换函数足够高效。标签编码保存图像时最好将类别标签也单独保存为一个.npy文件或csv文件与图像路径对应。更常见的做法是使用深度学习框架如PyTorch的ImageFolder能识别的目录结构即每个类别的图像放在一个子文件夹下框架会自动推断标签。GASF vs GADF没有绝对的好坏。我的经验是对于冲击特征明显的故障如点蚀GADF有时能提供更清晰的边缘信息。可以尝试将两者作为两个通道构建一个“2通道图像”输入CNN或者分别训练模型然后集成。6. 常见问题与排查技巧实录在实际复现和实验过程中你几乎一定会遇到下面这些问题。这里记录了我的排查思路和解决方法。6.1 生成的GAF图像全是灰色没有纹理现象保存的PNG图像看起来是均匀的灰色或者纹理非常微弱。可能原因与排查数据标准化错误检查normalize_sample函数。如果输入信号本身非常平缓方差很小Z-score标准化后数值范围可能仍然很窄比如在[-0.1, 0.1]导致arccos计算后角度变化不大生成的GAF矩阵元素值非常接近。解决打印几个样本标准化后的最大值和最小值。也可以尝试改用最小-最大归一化到[-1,1]看看是否改善。GAF值缩放错误在保存为图像前需要将GAF矩阵的值线性映射到[0, 255]。如果映射前的GAF矩阵动态范围很小映射后就会集中在中灰值附近。解决打印GAF矩阵的min()和max()。对于GASF理论范围是[-1,1]但实际可能集中在某个小区间。可以尝试使用对比度拉伸例如用(gaf_matrix - gaf_matrix.min()) / (gaf_matrix.max() - gaf_matrix.min()) * 255。下采样过于激进如果downsampled_len设得太小原始信号的特征模式在下采样过程中被严重平滑掉了。解决增大下采样长度或者先不做下采样生成大图看看是否有纹理。6.2 模型训练准确率始终在“瞎猜”水平比如四分类准确率25%左右现象CNN模型训练损失不下降验证准确率徘徊在随机猜测水平。可能原因与排查标签错乱这是最致命也最隐蔽的错误。检查数据加载环节确保图像路径和标签的对应关系绝对正确。一个快速验证的方法是取出训练集中每个类别的几张图像显示出来并打印其标签人工判断是否匹配。数据泄露确保在切片生成样本时没有使用跨越不同.mat文件的连续索引。例如文件A的末尾和文件B的开头被切到了一个样本里这个样本的标签就无法定义。确保每个样本完全来源于同一个数据文件。训练/测试集划分不当如果按文件顺序划分可能导致训练集和测试集来自完全不同的工况比如训练集全是0HP负载测试集全是3HP负载模型无法泛化。必须在样本级别进行随机打乱和分层划分确保每个集合中各类别、各工况的数据比例大致相同。GAF转换函数有数值错误使用一个小型人造信号测试你的gramian_angular_field函数。例如输入一个简单的正弦波观察生成的GAF图像是否具有预期的周期性结构。与公开的GAF实现如pyts库的结果进行对比。6.3 处理速度太慢特别是实时生成GAF时现象数据预处理或训练时的数据加载成为瓶颈。优化技巧向量化操作确保gramian_angular_field函数中使用了NumPy的向量化操作如np.outer避免Python层面的for循环。预处理与缓存如非必要不要在每个epoch实时生成GAF。花费一些时间和磁盘空间进行预处理是值得的。使用更高效的数据加载如果使用PyTorch将预处理好的图像数据集放在SSD硬盘上并使用torch.utils.data.DataLoader的多个工作进程num_workers和内存固定pin_memoryTrue来加速数据从磁盘到GPU的传输。降低图像分辨率这是最有效的提速方法。尝试128x128甚至64x64的图像很多情况下对分类精度影响不大但计算量和内存占用呈平方级下降。6.4 不同故障类型的GAF图像看起来区别不大现象人眼难以区分正常和故障的GAF图像。可能原因与对策故障特征微弱早期微弱故障或故障尺寸很小时振动信号中的冲击成分可能被强烈的背景噪声淹没。GAF转换并不能创造信息它只是换了一种表示方式。对策考虑在GAF转换前先对信号进行降噪处理例如使用小波阈值去噪、自适应滤波等。需要更强大的特征提取器人眼难以区分不代表CNN学不到特征。可以继续训练模型并利用Grad-CAM等可视化工具查看CNN到底关注图像的哪些区域。如果CNN能学到区分性特征那就没问题。尝试其他时频图像表示如果GAF效果确实不佳可以对比其他方法如连续小波变换CWT生成的时频谱图、短时傅里叶变换STFT谱图、马尔可夫变迁场MTF等。有时不同的数据适合不同的表示方法。轴承故障诊断是一个理论与实践紧密结合的领域。读懂东南大学的这份代码不仅仅是理解几行Python更是要理解其背后将时序问题转化为视觉问题的思想以及数据预处理中每一步的工程考量。从数据集的解读、清洗、切片到GAF图像的生成和优化每一步都藏着影响最终模型效果的细节。希望这份超详细的解读和实录的经验能帮你避开我当年踩过的那些坑更顺畅地踏上基于GAF和深度学习的智能诊断研究之路。在实际项目中多可视化、多对比、从小规模实验开始迭代是最高效的策略。