python基于三階貝塞爾曲線的數(shù)據(jù)平滑算法-創(chuàng)新互聯(lián)

前言

成都創(chuàng)新互聯(lián)于2013年成立,先為寧明等服務(wù)建站,寧明等地企業(yè),進(jìn)行企業(yè)商務(wù)咨詢(xún)服務(wù)。為寧明企業(yè)網(wǎng)站制作PC+手機(jī)+微官網(wǎng)三網(wǎng)同步一站式服務(wù)解決您的所有建站問(wèn)題。

很多文章在談及曲線平滑的時(shí)候,習(xí)慣使用擬合的概念,我認(rèn)為這是不恰當(dāng)?shù)摹F交蟮那€,一定經(jīng)過(guò)原始的數(shù)據(jù)點(diǎn),而擬合曲線,則不一定要經(jīng)過(guò)原始數(shù)據(jù)點(diǎn)。

一般而言,需要平滑的數(shù)據(jù)分為兩種:時(shí)間序列的單值數(shù)據(jù)、時(shí)間序列的二維數(shù)據(jù)。對(duì)于前者,并非一定要用貝塞爾算法,僅用樣條插值就可以輕松實(shí)現(xiàn)平滑;而對(duì)于后者,不管是 numpy 還是 scipy 提供的那些插值算法,就都不適用了。

本文基于三階貝塞爾曲線,實(shí)現(xiàn)了時(shí)間序列的單值數(shù)據(jù)和時(shí)間序列的二維數(shù)據(jù)的平滑算法,可滿(mǎn)足大多數(shù)的平滑需求。

貝塞爾曲線

關(guān)于貝塞爾曲線的數(shù)學(xué)原理,這里就不討論了,直接貼出結(jié)論:

一階貝塞爾曲線



python基于三階貝塞爾曲線的數(shù)據(jù)平滑算法

二階貝塞爾曲線



python基于三階貝塞爾曲線的數(shù)據(jù)平滑算法

三階貝塞爾曲線



python基于三階貝塞爾曲線的數(shù)據(jù)平滑算法

算法描述

如果我們把三階貝塞爾曲線的 P0 和 P3 視為原始數(shù)據(jù),只要找到 P1 和 P2 兩個(gè)點(diǎn)(我們稱(chēng)其為控制點(diǎn)),就可以根據(jù)三階貝塞爾曲線公式,計(jì)算出 P0 和 P3 之間平滑曲線上的任意點(diǎn)。


現(xiàn)在,平滑問(wèn)題變成了如何計(jì)算兩個(gè)原始數(shù)據(jù)點(diǎn)之間的控制點(diǎn)的問(wèn)題。步驟如下:

第1步:綠色直線連接相鄰的原始數(shù)據(jù)點(diǎn),計(jì)算出個(gè)線段的中點(diǎn),紅色直線連接相鄰的中點(diǎn)


第2步:根據(jù)相鄰兩條綠色直線長(zhǎng)度之比,分割其中點(diǎn)之間紅色連線,標(biāo)記分割點(diǎn)


第3步:平移紅色連線,使其分割點(diǎn)與相對(duì)的原始數(shù)據(jù)點(diǎn)重合


第4步:調(diào)整平移后紅色連線的端點(diǎn)與原始數(shù)據(jù)點(diǎn)的距離,通??s減40%-80%


算法實(shí)現(xiàn)

# -*- coding: utf-8 -*-

import numpy as np

def bezier_curve(p0, p1, p2, p3, inserted):
 """
 三階貝塞爾曲線
 
 p0, p1, p2, p3 - 點(diǎn)坐標(biāo),tuple、list或numpy.ndarray類(lèi)型
 inserted  - p0和p3之間插值的數(shù)量
 """
 
 assert isinstance(p0, (tuple, list, np.ndarray)), u'點(diǎn)坐標(biāo)不是期望的元組、列表或numpy數(shù)組類(lèi)型'
 assert isinstance(p0, (tuple, list, np.ndarray)), u'點(diǎn)坐標(biāo)不是期望的元組、列表或numpy數(shù)組類(lèi)型'
 assert isinstance(p0, (tuple, list, np.ndarray)), u'點(diǎn)坐標(biāo)不是期望的元組、列表或numpy數(shù)組類(lèi)型'
 assert isinstance(p0, (tuple, list, np.ndarray)), u'點(diǎn)坐標(biāo)不是期望的元組、列表或numpy數(shù)組類(lèi)型'
 
 if isinstance(p0, (tuple, list)):
  p0 = np.array(p0)
 if isinstance(p1, (tuple, list)):
  p1 = np.array(p1)
 if isinstance(p2, (tuple, list)):
  p2 = np.array(p2)
 if isinstance(p3, (tuple, list)):
  p3 = np.array(p3)
 
 points = list()
 for t in np.linspace(0, 1, inserted+2):
  points.append(p0*np.power((1-t),3) + 3*p1*t*np.power((1-t),2) + 3*p2*(1-t)*np.power(t,2) + p3*np.power(t,3))
 
 return np.vstack(points)


def smoothing_base_bezier(date_x, date_y, k=0.5, inserted=10, closed=False):
 """
 基于三階貝塞爾曲線的數(shù)據(jù)平滑算法
 
 date_x  - x維度數(shù)據(jù)集,list或numpy.ndarray類(lèi)型
 date_y  - y維度數(shù)據(jù)集,list或numpy.ndarray類(lèi)型
 k   - 調(diào)整平滑曲線形狀的因子,取值一般在0.2~0.6之間。默認(rèn)值為0.5
 inserted - 兩個(gè)原始數(shù)據(jù)點(diǎn)之間插值的數(shù)量。默認(rèn)值為10
 closed  - 曲線是否封閉,如是,則首尾相連。默認(rèn)曲線不封閉
 """
 
 assert isinstance(date_x, (list, np.ndarray)), u'x數(shù)據(jù)集不是期望的列表或numpy數(shù)組類(lèi)型'
 assert isinstance(date_y, (list, np.ndarray)), u'y數(shù)據(jù)集不是期望的列表或numpy數(shù)組類(lèi)型'
 
 if isinstance(date_x, list) and isinstance(date_y, list):
  assert len(date_x)==len(date_y), u'x數(shù)據(jù)集和y數(shù)據(jù)集長(zhǎng)度不匹配'
  date_x = np.array(date_x)
  date_y = np.array(date_y)
 elif isinstance(date_x, np.ndarray) and isinstance(date_y, np.ndarray):
  assert date_x.shape==date_y.shape, u'x數(shù)據(jù)集和y數(shù)據(jù)集長(zhǎng)度不匹配'
 else:
  raise Exception(u'x數(shù)據(jù)集或y數(shù)據(jù)集類(lèi)型錯(cuò)誤')
 
 # 第1步:生成原始數(shù)據(jù)折線中點(diǎn)集
 mid_points = list()
 for i in range(1, date_x.shape[0]):
  mid_points.append({
   'start': (date_x[i-1], date_y[i-1]),
   'end':  (date_x[i], date_y[i]),
   'mid':  ((date_x[i]+date_x[i-1])/2.0, (date_y[i]+date_y[i-1])/2.0)
  })
 
 if closed:
  mid_points.append({
   'start': (date_x[-1], date_y[-1]),
   'end':  (date_x[0], date_y[0]),
   'mid':  ((date_x[0]+date_x[-1])/2.0, (date_y[0]+date_y[-1])/2.0)
  })
 
 # 第2步:找出中點(diǎn)連線及其分割點(diǎn)
 split_points = list()
 for i in range(len(mid_points)):
  if i < (len(mid_points)-1):
   j = i+1
  elif closed:
   j = 0
  else:
   continue
  
  x00, y00 = mid_points[i]['start']
  x01, y01 = mid_points[i]['end']
  x10, y10 = mid_points[j]['start']
  x11, y11 = mid_points[j]['end']
  d0 = np.sqrt(np.power((x00-x01), 2) + np.power((y00-y01), 2))
  d1 = np.sqrt(np.power((x10-x11), 2) + np.power((y10-y11), 2))
  k_split = 1.0*d0/(d0+d1)
  
  mx0, my0 = mid_points[i]['mid']
  mx1, my1 = mid_points[j]['mid']
  
  split_points.append({
   'start': (mx0, my0),
   'end':  (mx1, my1),
   'split': (mx0+(mx1-mx0)*k_split, my0+(my1-my0)*k_split)
  })
 
 # 第3步:平移中點(diǎn)連線,調(diào)整端點(diǎn),生成控制點(diǎn)
 crt_points = list()
 for i in range(len(split_points)):
  vx, vy = mid_points[i]['end'] # 當(dāng)前頂點(diǎn)的坐標(biāo)
  dx = vx - split_points[i]['split'][0] # 平移線段x偏移量
  dy = vy - split_points[i]['split'][1] # 平移線段y偏移量
  
  sx, sy = split_points[i]['start'][0]+dx, split_points[i]['start'][1]+dy # 平移后線段起點(diǎn)坐標(biāo)
  ex, ey = split_points[i]['end'][0]+dx, split_points[i]['end'][1]+dy # 平移后線段終點(diǎn)坐標(biāo)
  
  cp0 = sx+(vx-sx)*k, sy+(vy-sy)*k # 控制點(diǎn)坐標(biāo)
  cp1 = ex+(vx-ex)*k, ey+(vy-ey)*k # 控制點(diǎn)坐標(biāo)
  
  if crt_points:
   crt_points[-1].insert(2, cp0)
  else:
   crt_points.append([mid_points[0]['start'], cp0, mid_points[0]['end']])
  
  if closed:
   if i < (len(mid_points)-1):
    crt_points.append([mid_points[i+1]['start'], cp1, mid_points[i+1]['end']])
   else:
    crt_points[0].insert(1, cp1)
  else:
   if i < (len(mid_points)-2):
    crt_points.append([mid_points[i+1]['start'], cp1, mid_points[i+1]['end']])
   else:
    crt_points.append([mid_points[i+1]['start'], cp1, mid_points[i+1]['end'], mid_points[i+1]['end']])
    crt_points[0].insert(1, mid_points[0]['start'])
 
 # 第4步:應(yīng)用貝塞爾曲線方程插值
 out = list()
 for item in crt_points:
  group = bezier_curve(item[0], item[1], item[2], item[3], inserted)
  out.append(group[:-1])
 
 out.append(group[-1:])
 out = np.vstack(out)
 
 return out.T[0], out.T[1]


if __name__ == '__main__':
 import matplotlib.pyplot as plt
 
 x = np.array([2,4,4,3,2])
 y = np.array([2,2,4,3,4])
	
	plt.plot(x, y, 'ro')
 x_curve, y_curve = smoothing_base_bezier(x, y, k=0.3, closed=True)
 plt.plot(x_curve, y_curve, label='$k=0.3$')
 x_curve, y_curve = smoothing_base_bezier(x, y, k=0.4, closed=True)
 plt.plot(x_curve, y_curve, label='$k=0.4$')
 x_curve, y_curve = smoothing_base_bezier(x, y, k=0.5, closed=True)
 plt.plot(x_curve, y_curve, label='$k=0.5$')
 x_curve, y_curve = smoothing_base_bezier(x, y, k=0.6, closed=True)
 plt.plot(x_curve, y_curve, label='$k=0.6$')
 plt.legend(loc='best')
 
 plt.show()

分享文章:python基于三階貝塞爾曲線的數(shù)據(jù)平滑算法-創(chuàng)新互聯(lián)
標(biāo)題路徑:http://muchs.cn/article16/higdg.html

成都網(wǎng)站建設(shè)公司_創(chuàng)新互聯(lián),為您提供全網(wǎng)營(yíng)銷(xiāo)推廣、云服務(wù)器網(wǎng)站策劃、網(wǎng)站設(shè)計(jì)、網(wǎng)頁(yè)設(shè)計(jì)公司、軟件開(kāi)發(fā)

廣告

聲明:本網(wǎng)站發(fā)布的內(nèi)容(圖片、視頻和文字)以用戶(hù)投稿、用戶(hù)轉(zhuǎn)載內(nèi)容為主,如果涉及侵權(quán)請(qǐng)盡快告知,我們將會(huì)在第一時(shí)間刪除。文章觀點(diǎn)不代表本網(wǎng)站立場(chǎng),如需處理請(qǐng)聯(lián)系客服。電話:028-86922220;郵箱:631063699@qq.com。內(nèi)容未經(jīng)允許不得轉(zhuǎn)載,或轉(zhuǎn)載時(shí)需注明來(lái)源: 創(chuàng)新互聯(lián)

微信小程序開(kāi)發(fā)