Savitzky-Golay 滤波是一种基于局部多项式拟合的数据平滑技术，广泛用于信号处理和时间序列分析
·
Savitzky-Golay 滤波是一种基于局部多项式拟合的数据平滑技术，广泛用于信号处理和时间序列分析。它通过在滑动窗口内对数据进行多项式拟合，平滑噪声同时尽可能保留数据的形状（如峰值、拐点），非常适合你的示波器多通道数据（SetMultiChannelViewData 方法中的 List<double> 数据）处理需求。以下是 Savitzky-Golay 滤波的详细原理、实现方法、在你的代码中的集成方式，以及优化建议，附带中文解释和代码示例。
1. Savitzky-Golay 滤波原理核心思想：
- 在一个固定大小的滑动窗口内（通常为奇数点），对数据点拟合一个低阶多项式（通常为 2 或 3 阶）。
- 使用拟合多项式的值（通常是窗口中心点的值）替换原始数据点。
- 通过滑动窗口逐点处理整个数据集，实现平滑。
关键参数：
- 窗口大小 (
2m+
)：控制平滑程度，越大越平滑，但可能丢失细节。1 - 多项式阶数 (( n ))：控制拟合的灵活性，通常为 2 或 3，过高可能过拟合噪声。
优点：
- 保留数据的局部特征（如峰值、拐点），优于简单移动平均。
- 对噪声鲁棒，适合示波器数据（如漏电流）。
- 卷积实现高效，适合大数据量。
缺点：
- 计算复杂度高于移动平均（需要矩阵运算或预计算系数）。
- 边界点处理需要特殊逻辑（可能引入误差）。
- 参数（窗口大小、阶数）需根据数据特性调优。
适用场景：
- 你的场景中，示波器数据（lst[id]）可能是时间序列，包含噪声但需要保留峰值或快速变化。Savitzky-Golay 滤波适合平滑漏电流（LEAKOSC）数据，同时保持关键特征。
2. 实现方法Savitzky-Golay 滤波的实现通常有两种方式：
- 直接拟合：在每个窗口内进行多项式拟合，灵活但计算量大。
- 卷积方法：预计算卷积系数，将平滑过程简化为卷积，效率更高。
由于你的应用需要处理多通道数据（SortedDictionary<int, List<double>> lst），推荐使用卷积方法，并借助库（如 MathNet.Numerics）简化实现。如果不使用库，也可以通过预计算系数实现。2.1 使用 MathNet.Numerics 库MathNet.Numerics 是一个强大的数值计算库，支持 Savitzky-Golay 滤波。以下是实现代码，需通过 NuGet 安装 MathNet.Numerics。安装：bash
dotnet add package MathNet.Numerics
代码：csharp
using MathNet.Numerics;
using MathNet.Numerics.Data.Text;
List<double> SmoothSavitzkyGolay(List<double> data, int windowSize, int polynomialOrder)
{
if (data == null || data.Count < windowSize || windowSize < polynomialOrder + 1 || windowSize % 2 == 0)
{
Debug.WriteLine("Invalid parameters for Savitzky-Golay filter");
return new List<double>(data);
}
// 转换为 MathNet 向量
var input = Vector<double>.Build.Dense(data.ToArray());
// 应用 Savitzky-Golay 滤波
// MathNet.Numerics 没有直接的 Savitzky-Golay 方法，但可通过卷积实现
// 这里使用自定义卷积系数（见下文）或外部库扩展
var coefficients = GenerateSavitzkyGolayCoefficients(windowSize, polynomialOrder);
var smoothed = new double[data.Count];
int halfWindow = windowSize / 2;
for (int i = 0; i < data.Count; i++)
{
double sum = 0;
for (int j = -halfWindow; j <= halfWindow; j++)
{
int index = i + j;
if (index >= 0 && index < data.Count)
{
sum += coefficients[j + halfWindow] * data[index];
}
}
smoothed[i] = sum;
}
return smoothed.ToList();
}
// 生成 Savitzky-Golay 卷积系数
double[] GenerateSavitzkyGolayCoefficients(int windowSize, int polynomialOrder)
{
int halfWindow = windowSize / 2;
var coefficients = new double[windowSize];
var x = Vector<double>.Build.Dense(windowSize, i => i - halfWindow);
var vandermonde = Matrix<double>.Build.Dense(windowSize, polynomialOrder + 1, (i, j) => Math.Pow(x[i], j));
var gram = vandermonde.TransposeThisAndMultiply(vandermonde);
var inverseGram = gram.Inverse();
var projection = inverseGram * vandermonde.Transpose();
var firstRow = projection.Row(0);
coefficients = firstRow.ToArray();
double sum = coefficients.Sum();
for (int i = 0; i < coefficients.Length; i++)
{
coefficients[i] /= sum; // 归一化
}
return coefficients;
}
说明：
- windowSize 必须为奇数（如 5、7、9），polynomialOrder 通常为 2 或 3。
- GenerateSavitzkyGolayCoefficients 计算卷积系数，基于最小二乘法。
- 对边界点，直接使用原始数据（可优化，见下文）。
- MathNet.Numerics 提供高效的矩阵运算，适合多通道数据。
2.2 手动实现（无外部库）如果不使用库，可以预计算常用配置的卷积系数，简化实现。以下是手动实现的示例（以窗口大小 5、多项式阶数 2 为例）。代码：csharp
List<double> SmoothSavitzkyGolay(List<double> data, int windowSize = 5, int polynomialOrder = 2)
{
if (data == null || data.Count < windowSize || windowSize % 2 == 0)
{
return new List<double>(data);
}
// 预定义卷积系数（窗口大小 5，阶数 2）
double[] coefficients = { -0.085714285714286, 0.342857142857143, 0.485714285714286, 0.342857142857143, -0.085714285714286 };
int halfWindow = windowSize / 2;
var smoothed = new List<double>();
for (int i = 0; i < data.Count; i++)
{
double sum = 0;
int count = 0;
for (int j = -halfWindow; j <= halfWindow; j++)
{
int index = i + j;
if (index >= 0 && index < data.Count)
{
sum += coefficients[j + halfWindow] * data[index];
count++;
}
}
smoothed.Add(count > 0 ? sum : data[i]); // 边界点用原始值
}
return smoothed;
}
说明：
- 使用预计算的系数（可从文献或工具生成），避免实时矩阵运算。
- 边界点直接复制原始值，简单但可能不平滑。
- 适用于窗口大小和阶数固定的场景，性能高。
3. 集成到 SetMultiChannelViewData以下是将 Savitzky-Golay 滤波（基于手动实现）集成到你的 SetMultiChannelViewData 方法的完整代码。使用 MathNet.Numerics 的实现可类似替换。代码：csharp
private void SetMultiChannelViewData(List<double> times_src, object obj2, Oscilloscope oscilloscope)
{
if (times_src == null || times_src.Count == 0 || obj2 == null)
{
Debug.WriteLine("Error: Invalid input data");
return;
}
string ChartType = oscilloscope.Name.Split('-')[1];
double MultiMaxY = double.MinValue;
double MultiMinY = double.MaxValue;
oscilloscope.LoadingData = true;
List<double> times = new List<double>(times_src);
SortedDictionary<int, List<double>> ls = obj2 as SortedDictionary<int, List<double>>;
SortedDictionary<int, List<double>> lst = new SortedDictionary<int, List<double>>();
// Savitzky-Golay 参数
int windowSize = 7; // 窗口大小，建议奇数
int polynomialOrder = 2; // 多项式阶数
// 数据复制和平滑
foreach (int id in ls.Keys)
{
if (PosIndexList.Contains(id))
{
lst.Add(id, SmoothSavitzkyGolay(ls[id], windowSize, polynomialOrder));
}
}
if (oscilloscope.PosNameRelation != SelectControlTestUnitSelecter.TestChannel.Bib.BIBCoords)
{
oscilloscope.PosNameRelation = SelectControlTestUnitSelecter.TestChannel.Bib.BIBCoords;
}
if (oscilloscope.Channels.Count == 0)
{
SetReatangleProps(SelectControlTestUnitSelecter.TestChannel.Id);
SetOscopeChannels(ChartType);
}
string unit = DataTimeConvertor.AutoConvertByTime(times[times.Count - 1], out double t_unit);
string yUnit = "";
double yScale = 1;
if (oscilloscope.Channels.Count > 0)
{
foreach (Channel c in oscilloscope.Channels)
{
int J = int.Parse(c.ChannelProps.Legend);
if (lst.ContainsKey(J) && lst[J].Count > 0)
{
double CacheMaxY = lst[J].Max();
double CacheMinY = lst[J].Min();
MultiMaxY = Math.Max(MultiMaxY, CacheMaxY);
MultiMinY = Math.Min(MultiMinY, CacheMinY);
}
}
if (MultiMaxY == MultiMinY)
{
MultiMaxY = MultiMinY + 0.1 * Math.Abs(MultiMinY);
}
// 动态范围计算
double range = Math.Abs(MultiMaxY - MultiMinY);
double padding = range == 0 ? Math.Abs(MultiMaxY) * 0.1 : range * 0.2;
double yMax = MultiMaxY + padding;
double yMin = MultiMinY - padding;
if (yMax < yMin)
{
double temp = yMax;
yMax = yMin;
yMin = temp;
}
if (ChartType == "LEAKOSC")
{
yUnit = DataTimeConvertor.AutoConvertByIc(MultiMaxY, out double maxYvalue, out yScale);
oscilloscope.FE = RuntimeConfiguration.ICConvert == "Scientific";
}
oscilloscope.BeginUpdate();
try
{
foreach (Channel c in oscilloscope.Channels)
{
int J = int.Parse(c.ChannelProps.Legend);
if (lst.ContainsKey(J) && lst[J].Count != 0)
{
c.Clear();
c.ChangeAxisXUnit(t_unit);
if (ChartType == "LEAKOSC")
{
c.ChangeAxisYUnit(yScale);
}
c.Add(times, lst[J]);
c.SetXAxisLimit(times[times.Count - 1], times[0], 0);
c.SetYAxisLimit(yMax, yMin, 0);
}
}
}
finally
{
oscilloscope.EndUpdate();
}
}
if (ChartType == "LEAKOSC")
{
oscilloscope.YLabel = $"Leakage [{yUnit}]";
}
oscilloscope.XLabel = $"Time [{unit}]";
}
// Savitzky-Golay 滤波实现
List<double> SmoothSavitzkyGolay(List<double> data, int windowSize = 7, int polynomialOrder = 2)
{
if (data == null || data.Count < windowSize || windowSize % 2 == 0)
{
return new List<double>(data);
}
// 预定义系数（窗口大小 7，阶数 2）
double[] coefficients = { -0.095238095238095, 0.142857142857143, 0.285714285714286, 0.333333333333333, 0.285714285714286, 0.142857142857143, -0.095238095238095 };
int halfWindow = windowSize / 2;
var smoothed = new List<double>();
for (int i = 0; i < data.Count; i++)
{
double sum = 0;
int count = 0;
for (int j = -halfWindow; j <= halfWindow; j++)
{
int index = i + j;
if (index >= 0 && index < data.Count)
{
sum += coefficients[j + halfWindow] * data[index];
count++;
}
}
smoothed.Add(count > 0 ? sum : data[i]);
}
return smoothed;
}
说明：
- 使用窗口大小 7、阶数 2 的预计算系数，适合漏电流数据（微小值，需保留峰值）。
- 平滑后的数据存储在 lst[id]，直接用于图表渲染。
- 边界点使用原始值，简单但可优化（见下文）。
4. 优化建议为进一步提升 Savitzky-Golay 滤波在你的示波器数据处理中的效果，以下是优化建议：4.1 边界点处理问题：当前实现对边界点（前 halfWindow 和后 halfWindow 个点）直接使用原始值，可能导致平滑不连续。优化：
- 外推法：在边界外补充数据（如复制边界值或镜像数据），然后应用滤波。
- 缩小窗口：在边界处动态减小窗口大小。
代码示例（外推法）：csharp
List<double> SmoothSavitzkyGolay(List<double> data, int windowSize = 7, int polynomialOrder = 2)
{
if (data == null || data.Count < windowSize || windowSize % 2 == 0)
{
return new List<double>(data);
}
double[] coefficients = { -0.095238095238095, 0.142857142857143, 0.285714285714286, 0.333333333333333, 0.285714285714286, 0.142857142857143, -0.095238095238095 };
int halfWindow = windowSize / 2;
var smoothed = new List<double>();
// 外推数据
var extendedData = new List<double>();
for (int i = 0; i < halfWindow; i++) extendedData.Add(data[0]); // 左边界补值
extendedData.AddRange(data);
for (int i = 0; i < halfWindow; i++) extendedData.Add(data[data.Count - 1]); // 右边界补值
for (int i = 0; i < data.Count; i++)
{
double sum = 0;
for (int j = -halfWindow; j <= halfWindow; j++)
{
sum += coefficients[j + halfWindow] * extendedData[i + j + halfWindow];
}
smoothed.Add(sum);
}
return smoothed;
}
说明：
- 在数据两端补充边界值，确保所有点都能应用完整窗口。
- 提高边界平滑效果，曲线更连续。
4.2 参数调优问题：窗口大小和阶数需根据数据特性调整，固定值可能不适合所有场景。优化：
- 动态窗口：根据数据点密度或噪声水平调整窗口大小（如数据点多时增大窗口）。
- 自适应阶数：根据数据复杂性选择阶数（简单趋势用 2 阶，复杂用 3 阶）。
- 用户可调：允许用户通过 UI 设置参数。
建议：
- 漏电流数据（微小值，-0.000116 到 -0.001022665）：推荐窗口大小 5-11，阶数 2。
- 测试不同参数，观察图表是否保留峰值且平滑噪声。
4.3 性能优化问题：对于大数据量（>10万点）或多通道，滤波可能耗时。优化：
- 预计算系数：如当前实现，使用固定系数避免实时矩阵运算。
- 下采样：在平滑前减少数据点数（参考上一回答）。
- 并行处理：对多通道数据并行应用滤波。
代码示例（下采样结合平滑）：csharp
List<double> DownsampleData(List<double> data, int targetPoints)
{
if (data.Count <= targetPoints) return new List<double>(data);
List<double> sampled = new List<double>();
double step = (double)data.Count / targetPoints;
for (int i = 0; i < targetPoints; i++)
{
sampled.Add(data[(int)(i * step)]);
}
return sampled;
}
// 在 SetMultiChannelViewData 中
int targetPoints = 1000; // 根据图表宽度
foreach (int id in ls.Keys)
{
if (PosIndexList.Contains(id))
{
var sampled = DownsampleData(ls[id], targetPoints);
lst.Add(id, SmoothSavitzkyGolay(sampled, windowSize, polynomialOrder));
}
}
times = DownsampleData(times, targetPoints);
说明：
- 下采样减少数据点数，降低滤波计算量。
- 确保采样后仍保留主要趋势。
4.4 异常值处理问题：漏电流数据可能包含异常值（如尖峰），影响平滑效果。优化：
- 预处理：在平滑前剔除极端值（如使用分位数）。
- 中值滤波：先应用中值滤波去除尖峰，再用 Savitzky-Golay 平滑。
代码示例（中值滤波预处理）：csharp
List<double> MedianFilter(List<double> data, int windowSize = 3)
{
if (data == null || data.Count < windowSize) return new List<double>(data);
var filtered = new List<double>();
int halfWindow = windowSize / 2;
for (int i = 0; i < data.Count; i++)
{
var window = new List<double>();
for (int j = -halfWindow; j <= halfWindow; j++)
{
int index = i + j;
if (index >= 0 && index < data.Count) window.Add(data[index]);
}
filtered.Add(window.OrderBy(x => x).ElementAt(window.Count / 2));
}
return filtered;
}
// 在 SmoothSavitzkyGolay 前
lst[id] = SmoothSavitzkyGolay(MedianFilter(ls[id], 3), windowSize, polynomialOrder);
说明：
- 中值滤波有效去除尖峰，保护 Savitzky-Golay 的平滑效果。
5. 测试与调优为确保 Savitzky-Golay 滤波在你的示波器图表中效果最佳，建议以下测试：
- 参数测试：
- 窗口大小：测试 5、7、9、11，观察平滑程度和细节保留。
- 阶数：测试 2 和 3，确认是否保留峰值。
- 示例：对数据 [-0.001022665, -0.000116, ...] 测试，检查曲线是否平滑且峰值清晰。
- 边界效果：
- 检查图表边界处是否平滑连续，使用外推法优化。
- 性能测试：
- 测试大数据量（>10万点）和多通道（>10个）场景，记录滤波时间。
- 结合下采样，目标保持滤波时间 <100ms。
- 视觉效果：
- 比较原始数据和平滑后数据的图表，确保噪声减少且趋势清晰。
- 验证多通道曲线是否一致（使用相同参数）。
6. 总结
- 优势：Savitzky-Golay 滤波适合你的漏电流数据，能平滑噪声同时保留峰值和快速变化，优于简单移动平均。
- 实现：
- 推荐使用 MathNet.Numerics 实现通用滤波，或手动实现固定参数以提高性能。
- 集成到 SetMultiChannelViewData，对 lst[id] 应用滤波。
- 优化：
- 边界处理：使用外推法。
- 性能：结合下采样和预计算系数。
- 异常值：预处理尖峰。
- 参数建议：
- 窗口大小：7-11
- 阶数：2
- 动态调整根据数据特性。
通过这些实现和优化，你的示波器图表将显示更平滑、清晰的曲线，提升数据分析效果。如果需要更详细的库配置、系数计算或其他平滑技术比较，请提供更多细节！
DAMO开发者矩阵，由阿里巴巴达摩院和中国互联网协会联合发起，致力于探讨最前沿的技术趋势与应用成果，搭建高质量的交流与分享平台，推动技术创新与产业应用链接，围绕“人工智能与新型计算”构建开放共享的开发者生态。
更多推荐
12
0- 0
已为社区贡献123条内容
所有评论(0)