用N点线性平滑法绘制曲线 木灵的炼金工作室

为了应对开学之后的实验报告,解放双手和脑袋去打代码,最近在开发面向不怎么喜欢的材料专业课的实验的实验数据处理工具。这些实验最喜欢玩儿的就是给几个非常离散的实验数据点,让你“作平滑曲线”。听起来很简单,平常手绘的时候没有很注意这个具体的绘制方式,但是现在要写程序,就需要一个算法上的描述。

线性插值

考虑我们的输入只有几个很离散的点,所以先补全其中间的点:

class smooth {
private:
	vector<pair<double, double>> origin;
	unsigned int resolution; //分辨率,每俩点之间插几个
	unsigned int smoothLevel; //线性平滑取样数
	unsigned int smoothTimes; //线性平滑次数

	vector<double> resX;
	vector<double> resY;

	void linearInsert();
	void linearSmooth();
public:
	explicit smooth(vector<pair<double, double>> &ori, unsigned int _reso = 25, unsigned int smoothLevel = 3, unsigned int smoothTime = 20);
	pair<vector<double>, vector<double>> runSmooth();
};

void smooth::linearInsert() {
	if (origin.size() <= 1) {
		throw("1个点平滑**呢");
	}
	for (int i = 0; i < origin.size() - 1; i++) {
		for (int j = 0; j < resolution + 1; j++) {
			resX.push_back((origin[i + 1].first - origin[i].first) * ((double)j / (double)(resolution + 1)) + origin[i].first);
			resY.push_back((origin[i + 1].second - origin[i].second) * ((double)j / (double)(resolution + 1)) + origin[i].second);
		}
	}
}

很直观的线性插值算法,在每两个相邻的离散点之间补入resolution个直线上的点。
插值完成后再用我之前写的绘图工具的魔改版本绘制其图形:

//测试用例:
vector<pair<double, double>> tv = {
	{1, -5},
	{2, 5},
	{3, -5},
	{4, 5},
	{5, -5},
	{6, 5}
};

插值数量为25的线性插值绘图
线性插值1

一点都不光滑。所以使用一个简单的平滑算法来处理。

移动平均 / N点线性平滑

移动平均是一个很常用也很直观的图像平滑算法,它对于图像上尖刺的去除十分有效。
它的实现就是把数据中每一点的值替换为其周围值的平均:看代码

void smooth::linearSmooth() {
	vector<double> tempRes(resY);
	for (int i = 0; i < resY.size(); i++) {
		double res = 0;
		for (int j = i > smoothLevel ? i - smoothLevel : 0; j <= min(i + smoothLevel, resY.size() - 1); j++) {
			res += resY[j];
		}
		tempRes[i] = res / (double)min((smoothLevel * 2 + 1), min(i + smoothLevel + 1, resY.size() - i + smoothLevel));
	}
	resY.~vector();
	resY = tempRes;
}

对上述图形执行一次这样的算法,得到的图形是:
线性插值2
其峰值被显著平滑化。

执行次数越多,平滑化越显著,但同时信息丢失也越严重(由信息论的一些基础知识可以很直观地得出这一结论)。比如。随着执行次数的增大:

执行5次:
线性插值3

执行20次:
线性插值4

执行50次:
线性插值5

执行100次:
线性插值6

执行250次:
线性插值7

执行500次:
线性插值8

算法是线性时空的,但是随着执行次数的增多,图形的失真会越发严重。所以建议这样的算法执行1-10次即可。
这个算法对于平滑峰值有着很显著的作用:即使运行一次,图形中的连续不可导点就会被磨去。这个点两侧的导数差距越大,这个点越容易被去除。
但是它对于“轻微的”不可导点效果不佳,以及被线性平滑处理的曲线不可以精确地穿过给出的离散点,考虑未来使用贝塞尔曲线。

另外,可以使用二项式系数来为平均值赋权以遏制图像的失真。

我们生造一组实验数据(意味深):

vector<pair<double, double>> tv = {
	{0, 0},
	{1, 5},
	{2, 25},
	{3, 20},
	{4, 18},
	{5, 17},
	{6, 14},
	{7, 3},
	{8, -6},
	{9, -3},
	{10, 3},
	{11, 4},
	{12, 9}
};

0次线性平滑
线性插值9

10次线性平滑 线性插值10

效果还不错嘛。

最后附上核心的代码:

#pragma once
#include "message.h"
#include <iostream>
#include <vector>
#include <algorithm>
#include <utility>

using namespace std;

class smooth {
private:
	vector<pair<double, double>> origin;
	unsigned int resolution; //分辨率,每俩点之间插几个
	unsigned int smoothLevel; //线性平滑取样数
	unsigned int smoothTimes; //线性平滑次数

	vector<double> resX;
	vector<double> resY;

	void linearInsert();
	void linearSmooth();
public:
	explicit smooth(vector<pair<double, double>> &ori, unsigned int _reso = 25, unsigned int smoothLevel = 3, unsigned int smoothTime = 3);
	pair<vector<double>, vector<double>> runSmooth();
};

smooth::smooth(vector<pair<double, double>> &ori, unsigned int _reso, unsigned int _smoLe, unsigned int _smoTi) 
	: origin(ori), resolution(_reso), smoothLevel(_smoLe), smoothTimes(_smoTi) {
	sort(origin.begin(), origin.end(), [&](pair<double, double> a, pair<double, double> b) {
		return a.first < b.first;
	});
}

void smooth::linearInsert() {
	if (origin.size() <= 1) {
		throw("1个点平滑**呢");
	}
	for (int i = 0; i < origin.size() - 1; i++) {
		for (int j = 0; j < resolution + 1; j++) {
			resX.push_back((origin[i + 1].first - origin[i].first) * ((double)j / (double)(resolution + 1)) + origin[i].first);
			resY.push_back((origin[i + 1].second - origin[i].second) * ((double)j / (double)(resolution + 1)) + origin[i].second);
		}
	}
}

void smooth::linearSmooth() {
	vector<double> tempRes(resY);
	for (int i = 0; i < resY.size(); i++) {
		double res = 0;
		for (int j = i > smoothLevel ? i - smoothLevel : 0; j <= min(i + smoothLevel, resY.size() - 1); j++) {
			res += resY[j];
		}
		tempRes[i] = res / (double)min((smoothLevel * 2 + 1), min(i + smoothLevel + 1, resY.size() - i + smoothLevel));
	}
	resY.~vector();
	resY = tempRes;
}

pair<vector<double>, vector<double>> smooth::runSmooth() {
	this->linearInsert();
	for (int i = 0; i < smoothTimes; i++){
		this->linearSmooth();
	}
	return {resX, resY};
}


Copyright AmachiInori 2017-2021. All Right Reserved.
Powered By Jekyll.
amachi.com.cn