110 lines
3.8 KiB
Python
110 lines
3.8 KiB
Python
import numpy as np
|
||
import pandas as pd
|
||
from pathlib import Path
|
||
from scipy.interpolate import interp1d, CubicSpline, UnivariateSpline, PchipInterpolator, Akima1DInterpolator
|
||
from statsmodels.tsa.stattools import acf
|
||
from scipy.signal import find_peaks
|
||
import matplotlib.pyplot as plt
|
||
|
||
# 1. 读取 CSV 文件并提取指定列
|
||
def read_csv_column(file_path, column_name):
|
||
if isinstance(file_path, str):
|
||
file_path = Path(file_path)
|
||
df = pd.read_csv(file_path)
|
||
if column_name not in df.columns:
|
||
raise ValueError(f"Column '{column_name}' not found in CSV file.")
|
||
return df[column_name].values
|
||
|
||
# 2. 计算自相关函数（ACF）
|
||
def compute_acf(data):
|
||
n = len(data)
|
||
nlags = n // 2
|
||
acf_values = acf(data, nlags=nlags, fft=True)
|
||
return acf_values[:nlags] # 截断为原长度的一半
|
||
|
||
# 3. 插值函数
|
||
"""
|
||
- 'linear' : 线性插值
|
||
- 'cubic' : 三次样条插值 (CubicSpline)
|
||
- 'spline' : 样条插值 (UnivariateSpline)
|
||
- 'pchip' : 分段三次 Hermite 插值 (PchipInterpolator)
|
||
- 'akima' : Akima 插值 (Akima1DInterpolator)
|
||
- 'quadratic' : 二次插值
|
||
"""
|
||
def interpolate_data(data, method='cubic', num_points=100, **kwargs):
|
||
x = np.arange(len(data))
|
||
|
||
if method == 'linear':
|
||
interp_func = interp1d(x, data, kind='linear', fill_value="extrapolate", **kwargs)
|
||
elif method == 'cubic':
|
||
interp_func = CubicSpline(x, data, bc_type='natural', **kwargs)
|
||
elif method == 'spline':
|
||
# UnivariateSpline 默认使用三次样条，但可以指定 degree=s
|
||
interp_func = UnivariateSpline(x, data, s=kwargs.get('s', None), k=kwargs.get('k', 3), **kwargs)
|
||
elif method == 'pchip':
|
||
interp_func = PchipInterpolator(x, data, **kwargs)
|
||
elif method == 'akima':
|
||
interp_func = Akima1DInterpolator(x, data, **kwargs)
|
||
elif method == 'quadratic':
|
||
interp_func = interp1d(x, data, kind='quadratic', fill_value="extrapolate", **kwargs)
|
||
else:
|
||
raise ValueError(f"Unsupported interpolation method: {method}. "
|
||
f"Supported methods: linear, cubic, spline, pchip, akima, quadratic, quintic")
|
||
|
||
x_new = np.linspace(0, len(data) - 1, num_points)
|
||
y_new = interp_func(x_new)
|
||
return y_new
|
||
|
||
# 4. 查找极大值点
|
||
def find_maxima(data):
|
||
peaks, _ = find_peaks(data, height=np.mean(data))
|
||
return peaks.tolist(), data[peaks].tolist()
|
||
|
||
# 5. 主流程
|
||
if __name__ == "__main__":
|
||
# 示例文件路径和列名
|
||
file_path = Path("尺寸标定板/ver-1.csv") # 假设当前目录下有 data.csv
|
||
column_name = "Gray_Value"
|
||
|
||
# Step 1: 读取 CSV 文件
|
||
data = read_csv_column(file_path, column_name)
|
||
print(len(data))
|
||
|
||
|
||
# Step 2: 计算自相关函数
|
||
acf_data = compute_acf(data)
|
||
print(len(acf_data))
|
||
|
||
# Step 3: 插值自相关函数
|
||
interpolated_acf = interpolate_data(acf_data, method='cubic', num_points=len(acf_data) * 1000)
|
||
print(len(interpolated_acf))
|
||
|
||
# Step 4: 查找极大值
|
||
peaks_indices, peaks_values = find_maxima(interpolated_acf)
|
||
|
||
# Step 5: 可视化
|
||
plt.figure(figsize=(12, 10))
|
||
|
||
# 原始数据
|
||
plt.subplot(3, 1, 1)
|
||
plt.plot(data, label='Original Data')
|
||
plt.title('Original Data')
|
||
plt.grid(True)
|
||
|
||
# 自相关函数
|
||
plt.subplot(3, 1, 2)
|
||
plt.plot(acf_data, label='ACF', marker='o')
|
||
plt.title('Autocorrelation Function (Truncated)')
|
||
plt.grid(True)
|
||
|
||
# 插值后的自相关函数及极大值
|
||
plt.subplot(3, 1, 3)
|
||
x_interpolated = np.linspace(0, len(acf_data) - 1, len(interpolated_acf))
|
||
print(x_interpolated[peaks_indices])
|
||
plt.plot(x_interpolated, interpolated_acf, label='Interpolated ACF', color='tab:orange')
|
||
plt.plot(x_interpolated[peaks_indices], peaks_values, 'rx', markersize=10, label='Maxima')
|
||
plt.title('Interpolated ACF with Maxima')
|
||
plt.grid(True)
|
||
plt.tight_layout()
|
||
plt.legend()
|
||
plt.show() |