【问题标题】:Curvature of set of points点集的曲率
【发布时间】:2018-05-28 19:16:40
【问题描述】:

data file

import matplotlib.pylab as plt
import numpy as np

#initial data
data=np.loadtxt('profile_nonoisebigd02.txt')
x=data[:,0]
y=data[:,1]

#first derivatives 
dx= np.gradient(data[:,0])
dy = np.gradient(data[:,1])

#second derivatives 
d2x = np.gradient(dx)
d2y = np.gradient(dy)

#calculation of curvature from the typical formula
curvature = np.abs(dx * d2y - d2x * dy) / (dx * dx + dy * dy)**1.5

谁能帮我看看我的曲率哪里出了问题? 这组点给了我一条抛物线,但曲率不是我所期望的。

【问题讨论】:

  • 好像两张图是一样的。
  • 我的代码得到了预期的结果。你能贴出你用来创建这两个图的 sn-p 吗?
  • plt.scatter(x,y) 用于初始图,plt.scatter(x,curvature) 用于曲率图。
  • 仍然无法重现您的问题。考虑发布data 的子集,可能类似于data[::data.shape[0]//25]
  • @YakymPirozhenko 我已经添加了整个数据文件。 Btw 是 d2x, d2y 计算正确吗?

标签: python curve


【解决方案1】:

您的数据似乎不够平滑;我使用 pandas 通过滚动方式替换不同值窗口大小的 x、y、dx、dy、d2x、d2y 和曲率。随着窗口大小的增加,曲率开始看起来越来越像您期望看到的平滑抛物线(图例给出了窗口大小):

作为参考,这里是您的原始数据图:

用于创建平滑帧的代码:

def get_smooth(smoothing=10, return_df=False):
    data=np.loadtxt('profile_nonoisebigd02.txt')

    if return_df:
        return pd.DataFrame(data)

    df = pd.DataFrame(data).sort_values(by=0).reset_index(drop=True).rolling(smoothing).mean().dropna()

    # first derivatives
    df['dx'] = np.gradient(df[0])
    df['dy'] = np.gradient(df[1])

    df['dx'] = df.dx.rolling(smoothing, center=True).mean()
    df['dy'] = df.dy.rolling(smoothing, center=True).mean()

    # second derivatives
    df['d2x'] = np.gradient(df.dx)
    df['d2y'] = np.gradient(df.dy)

    df['d2x'] = df.d2x.rolling(smoothing, center=True).mean()
    df['d2y'] = df.d2y.rolling(smoothing, center=True).mean()


    # calculation of curvature from the typical formula
    df['curvature'] = df.eval('abs(dx * d2y - d2x * dy) / (dx * dx + dy * dy) ** 1.5')
    # mask = curvature < 100

    df['curvature'] = df.curvature.rolling(smoothing, center=True).mean()

    df.dropna(inplace=True)
    return df[0], df.curvature

【讨论】:

    猜你喜欢
    • 2013-06-19
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 2014-06-22
    相关资源
    最近更新 更多