2013-07-29 1 views
5

주기 인장 시험에서 데이터를 분석합니다. 입력으로 사용되는 거대한 x와 y 값의 목록입니다. 소재가 딱딱 해지거나 부드럽게되는지 설명하려면 각 사이클 루프의 파란색 슬로프를 가져와야합니다. xy 데이터 포인트 그래프와 numpy의 모든 교차점을 제거 하시겠습니까?

tensile_test

경사의 낮은 지점을 얻기

slope

아동 축제이지만, 상단 하나, 그것은 도전이다. 나는 각 루프의 로컬 최대 아래 몇 가지 포인트로 루프를 슬라이스하고 에서 빨간 선을, 지금까지이 방법을 만들었습니다

the_challenge

포인트의 수를 hardnumbered. 빨간색 선의 근접은 poly1d(polyfit(x1,x2,1))fsolve에 의해 만들어지고 교점을 얻는 데 사용됩니다. 그러나 점의 분포가 항상 같지 않기 때문에 항상 올바르게 작동하지는 않습니다.

문제가 올바르게 교차하는 선 간격 두 (빨강) 정의 방법이다. 위의 그림에서 평균 경사와 함께 3 번의 실험이 있습니다. 나는 최선의 방법이 아니라고 결정하면서 각 루프마다 4 개의 가장 가까운 점을 찾으려고 며칠을 보냈다. 그리고 마침내, 저는 여기 stackowerflow에서 끝났습니다.

원하는 출력은 교점의 대략적인 좌표를 가진 목록입니다. 재생하려는 경우 here은 곡선 (0, [[xvals], [yvals]])의 데이터입니다. Theese 쉽게 당신은 같은 X 값에서 모든 라인의 y 값을 리샘플링 scipy에서 interp1d 기능을 사용하여 매우 간단하게이 작업을 수행 할 수

import csv 
import sys 
csv. field_size_limit(sys.maxsize)  

csvfile = 'data.csv' 
tc_data = {} 
for key, val in csv.reader(open(csvfile, "r")): 
    tc_data[key] = val 
for key in tc_data: 
    tc = eval(tc_data[key]) 

x = tc[0] 
y = tc[1] 
+0

을 나를 위해 일하고있다. – gggg

+1

matplotlib로 확대 효과를 얻은 방법이 궁금합니다 –

답변

6

이 조금 과잉, 그러나 당신의 교차점을 찾을 수있는 적절한 방법이 될 수있다, 첫 번째 덩어리에서 어떤 세그먼트가 어떤 세그먼트와 교차하는지 확인하는 것입니다 두 번째 청크에서.

나 자신에게 약간의 데이터를 쉽게하는 prolate cycloid의 조각을 만들려고하고, 그리고 here 유사하게 감소로 증가하는 넘겼습니다 Y 좌표 곳에 장소를 찾을 것입니다 :

a, b = 1, 2 
phi = np.linspace(3, 10, 100) 
x = a*phi - b*np.sin(phi) 
y = a - b*np.cos(phi) 
y_growth_flips = np.where(np.diff(np.diff(y) > 0))[0] + 1 

plt.plot(x, y, 'rx') 
plt.plot(x[y_growth_flips], y[y_growth_flips], 'bo') 
plt.axis([2, 12, -1.5, 3.5]) 
plt.show() 

enter image description here

P0에서 P1까지가는 세그먼트와 Q0에서 Q1까지가는 세그먼트가있는 경우 벡터 수식을 해결하면 교차점을 찾을 수 있습니다3210이고, 두 세그먼트는 실제로 st이 모두 [0, 1] 인 경우 교차합니다.모든 세그먼트이 밖으로 시도 : 내 시스템에서

x_down = x[y_growth_flips[0]:y_growth_flips[1]+1] 
y_down = y[y_growth_flips[0]:y_growth_flips[1]+1] 
x_up = x[y_growth_flips[1]:y_growth_flips[2]+1] 
y_up = y[y_growth_flips[1]:y_growth_flips[2]+1] 

def find_intersect(x_down, y_down, x_up, y_up): 
    for j in xrange(len(x_down)-1): 
     p0 = np.array([x_down[j], y_down[j]]) 
     p1 = np.array([x_down[j+1], y_down[j+1]]) 
     for k in xrange(len(x_up)-1): 
      q0 = np.array([x_up[k], y_up[k]]) 
      q1 = np.array([x_up[k+1], y_up[k+1]]) 
      params = np.linalg.solve(np.column_stack((p1-p0, q0-q1)), 
            q0-p0) 
      if np.all((params >= 0) & (params <= 1)): 
       return p0 + params[0]*(p1 - p0) 

>>> find_intersect(x_down, y_down, x_up, y_up) 
array([ 6.28302264, 1.63658676]) 

crossing_point = find_intersect(x_down, y_down, x_up, y_up) 
plt.plot(crossing_point[0], crossing_point[1], 'ro') 
plt.show() 

enter image description here

, 이것은 아마도 충분한 약 20 초고속없는 초당 교차로,하지만를 처리 할 수있는 것은 모든 이제 다음 그래프를 분석합니다. 당신은 2 × 2 선형 시스템의 솔루션을 벡터화하여 일을 spped 할 수 있습니다 :

def find_intersect_vec(x_down, y_down, x_up, y_up): 
    p = np.column_stack((x_down, y_down)) 
    q = np.column_stack((x_up, y_up)) 
    p0, p1, q0, q1 = p[:-1], p[1:], q[:-1], q[1:] 
    rhs = q0 - p0[:, np.newaxis, :] 
    mat = np.empty((len(p0), len(q0), 2, 2)) 
    mat[..., 0] = (p1 - p0)[:, np.newaxis] 
    mat[..., 1] = q0 - q1 
    mat_inv = -mat.copy() 
    mat_inv[..., 0, 0] = mat[..., 1, 1] 
    mat_inv[..., 1, 1] = mat[..., 0, 0] 
    det = mat[..., 0, 0] * mat[..., 1, 1] - mat[..., 0, 1] * mat[..., 1, 0] 
    mat_inv /= det[..., np.newaxis, np.newaxis] 
    import numpy.core.umath_tests as ut 
    params = ut.matrix_multiply(mat_inv, rhs[..., np.newaxis]) 
    intersection = np.all((params >= 0) & (params <= 1), axis=(-1, -2)) 
    p0_s = params[intersection, 0, :] * mat[intersection, :, 0] 
    return p0_s + p0[np.where(intersection)[0]] 

네, 지저분한, 그러나 그것은 작동, 그래서 × 100 배 빠르게 수행합니다 당신의 연결의 어느

find_intersect(x_down, y_down, x_up, y_up) 
Out[67]: array([ 6.28302264, 1.63658676]) 

find_intersect_vec(x_down, y_down, x_up, y_up) 
Out[68]: array([[ 6.28302264, 1.63658676]]) 

%timeit find_intersect(x_down, y_down, x_up, y_up) 
10 loops, best of 3: 66.1 ms per loop 

%timeit find_intersect_vec(x_down, y_down, x_up, y_up) 
1000 loops, best of 3: 375 us per loop 
+0

제이미, 맥주 빚이 있습니다. 나에게 주소를 보내면 즉시 발송할 것입니다! – ptaeck

+0

내가 TypeError를 얻고 있습니다 : int는'find_intersect_vec'의 prolest와 함께'intersection = np.all ((params> = 0) & (params <= 1), axis = (- 1, -2))'라인에 필요합니다. 사이클로이드 예제 데이터. – ptaeck

+0

어떤 numpy 버전을 사용하고 있습니까? 다중 축, 1.7에서 도입 된 곳, 이전 버전을 사용한다면 두 개의 호출을'np.all'에 중첩시켜야 할 수도 있습니다.'intersection = np.all (np.all ((params> = 0) & (params <= 1), axis = -1), axis = -1)'1.6에서도 같은 작업을 수행해야합니다. – Jaime