2016-05-31 6 views
1

Я застреваю с запутанной проблемой. Вот немного фона:Почему мой полный круг не на 360 градусов?

Я работаю над qgis/python с координатами в Lambert93: одна центральная точка (моя диктофонная клавиша) и несколько других точек, тяготеющих вокруг нее. Для упрощения кода я поставил только один пример:

import numpy as np 
import math 

dict = {(355385,6.68906e+06): [(355277,6.68901e+06), (355501,6.68912e+06), (355364,6.6891e+06), (355277,6.68901e+06)]} 

for key, values in dict.iteritems(): 
    anglist =[] 
    print key 
    i=0 
    j=1 
    for sides in values[:-1]: 

     A = np.array(dict[key][i]) 
     B = np.array(key) 
     C = np.array(dict[key][j]) 

     BA = A - B 
     BC = C - B 

     cosine_angle = np.vdot(BA, BC)/(np.linalg.norm(BA) * np.linalg.norm(BC)) 
     angle = (np.degrees(np.arccos(cosine_angle))) 

     i+=1 
     j+=1 

     anglist.append(angle) 
     s = sum(anglist) 
    dict[key]= [values, anglist, s] 
print dict 

Результаты:

{(355385, 6689060.0): [[(355277, 6689010.0), (355501, 6689120.0), (355364, 6689100.0), (355277, 6689010.0)], [177.4925133253854, 90.349597027985112, 87.142916297400205], 354.98502665077069]} 

Как вы можете видеть, сумма = 354. У меня есть большой набор данных и иногда получить правильный 360, но по большей части я этого не делаю. Но во всей логике, обернув одну точку и завершив расчет, где она началась, единственный результат, который я должен получить, равен 360.

Я пробовал второй способ, чтобы посмотреть, не были ли и angle проблема:

from math import sqrt 
from math import acos 
import numpy 
def angle(a, b, c): 

    # Create vectors from points 
    ba = [ aa-bb for aa,bb in zip(a,b) ] 
    bc = [ cc-bb for cc,bb in zip(c,b) ] 

# Normalize vector 
    nba = sqrt (sum ((x**2.0 for x in ba))) 
    ba = [ x/nba for x in ba ] 

    nbc = sqrt (sum ((x**2.0 for x in bc))) 
    bc = [ x/nbc for x in bc ] 

# Calculate scalar from normalized vectors 
    scale = sum ((aa*bb for aa,bb in zip(ba,bc))) 

# calculate the angle in radian 
    angle = numpy.degrees(acos(scale)) 
    return angle 

print angle((355277,6.68901e+06),(355385,6.68906e+06), (355501,6.68912e+06)) 
print angle((355501,6.68912e+06),(355385,6.68906e+06), (355364,6.6891e+06)) 
print angle((355364,6.6891e+06),(355385,6.68906e+06), (355277,6.68901e+06)) 

Но результаты по-прежнему:

177.492513325 
90.349597028 
87.1429162974 

Так что я думаю, что мы можем пересечь математику из проблемы ... Так что одна возможность проблема с тем, как QGIS (или питона?) управляет координатами. Как я могу обойти это?

Я должен сказать, коды в основном так же, как here, here и here

+1

В последнем примере, вы, вероятно, есть один угол, который> 180 градусов. Проблема в том, что стандартная ветвь функции acos() 'может только возвращать значения из [0, pi/2]. Скорее всего, первый угол должен быть 182,51 градуса - тогда цифры суммируются хорошо. Но сейчас у меня нет времени проверять цифры. Чтобы исправить это, вам нужна дополнительная контрольная логика. –

ответ

0

Благодаря индикации Hellmar, вот это рабочий код. Я сделал две вещи: дал все центральные точки (0,0) в своем самолете, вычитал расстояние между центральной точкой и (0,0) от всех других связанных точек. Тогда я смог использовать atan2, потому что все углы, которые мне нужно вычислить, были в (0,0). Atan2 очень практичен в этом контексте, потому что для вычисления угла нужно всего две точки: измеряемый угол всегда равен 0,0. Это, вероятно, означает, что мне не нужно устанавливать мою центральную точку как 0,0, так как ее просто выгоняют из уравнения. Вот мой код. Любые дальнейшие советы приветствуются.

import numpy as np 
dictionary = {(355385,6.68906e+06): [(355277,6.68901e+06), (355501,6.68912e+06), (355364,6.6891e+06), (355277,6.68901e+06)], (355364,6.6891e+06): [(355261,6.68905e+06), (355385,6.68906e+06), (355481,6.68916e+06), (355340,6.68915e+06), (355261,6.68905e+06)], (355340,6.68915e+06): [(355238,6.68909e+06), (355364,6.6891e+06), (355452,6.68921e+06), (355238,6.68909e+06)]} 
def angle_between(p1, p2): 
    ang1 = np.arctan2(*p1[::-1]) 
    ang2 = np.arctan2(*p2[::-1]) 
    return np.rad2deg((ang1 - ang2) % (2 * np.pi)) 

zlist=[] 
newdict={} 
for key, values in dictionary.iteritems(): 
    xlist =[] 
    ylist = [] 
    i=0 
#print key 
#print key[0] 
    for sides in values: 

     A = dictionary[key][i][0] 
     B = key[0] 
     C = dictionary[key][i][1] 
     D= key[1] 
     E = (0.0, 0.0) 

     o1 = A-B 
     o2 = C-D 

     xlist.append(o1) 
     ylist.append(o2) 
     ziplist = zip(xlist, ylist) 
     i+=1 
    zlist.append(ziplist) 

#print dict[key][i][0] 

newdict=zlist 
print newdict 
angledict = [] 

for p in newdict: 
    i=0 
    j=1 
    print p 
    for q in p[:-1]: 
     A=p[i] 
     B=p[j] 
     print "A=",A 
     print "B=", B 

     angledict.append(angle_between(A,B)) 
     i+=1 
     j+=1 

print angledict 

Я взял angle_between функцию от here