Циркулярная интерполяция в Python

Question

У меня есть две системы, каждая из которых имеет датчик направления (0-360 градусов), но датчики могут давать совершенно разные значения в зависимости от ориентации каждой системы и линейности каждого датчика. У меня есть механическая ссылка, которую я могу использовать для создания таблицы, в которой каждая система фактически указывает. Это дает таблицу с тремя столбцами:

Physical SystemA SystemB 
-------- ------- ------- 
000.0  005.7 182.3 
005.0  009.8 178.4 
...  ...  ... 

С показаны только данные, мы можем видеть, что SYSTEMA не далеко от физической ссылки, но SYSTEMB составляет около 180 градусов прочь, и идет в обратном направлении (представьте, что он установлен вверх дном).

Мне нужно иметь возможность сопоставлять между всеми тремя значениями: если сообщение SystemA что-то находится в 105.7, мне нужно сообщить пользователю, какое именно физическое направление, а затем сообщить SystemB о том же местоположении. То же самое, если SystemB делает первоначальный отчет. И пользователь может запросить обе системы указать желаемое физическое направление, поэтому SystemA и SystemB необходимо будет указать, где указать.

Линейная интерполяция не сложна, но у меня возникают проблемы, когда данные идут в противоположных направлениях и являются модульными/циклическими.

Есть ли способ для Pythonic делать все эти сопоставления?

---

EDIT: Давайте сосредоточимся на наиболее сложном случае, когда мы имеем два спаренных списка значений:

A  B 
----- ----- 
    0.0 182.5 
10.0 172.3 
20.0 161.4 
...  ... 
170.0  9.7 
180.0 359.1 
190.0 348.2 
...  ... 
340.0 163.6 
350.0 171.8 

Скажем, списки поступают из двух различных радаров с указателями, которые не выровнены по Север или что-то еще, но мы вручную взяли вышеуказанные данные, перемещая цель вокруг и видя, где каждый радар должен был указать, чтобы увидеть это.

Когда радар А говорит: «У меня есть цель в 123.4!», Где мне нужно прицелиться в радар B, чтобы увидеть его? Если радар B находит цель, где я должен указывать радар A?

Список A Обертывание между последним и первым элементом, но список B располагается ближе к середине списка. Список A монотонно возрастает, а список B монотонно уменьшается. Обратите внимание, что размер степени на А, как правило, не совпадает с размером степени на В.

Есть простой интерполятор, который будет обернуть правильно, когда:

1. Интерполирующие из списка А в список B .

2. Интерполирующие из списка B в список A.

это нормально использовать два отдельных интерполятор инстанциацию, один для собирается в каждом направлении. Я предполагаю, что линейный (первый порядок) интерполятор в порядке, но я, возможно, захочу использовать более высокую или сплайнную интерполяцию в будущем.

В некоторых случаях тест:

• A = 356.7, В =?

• A = 179.2, B =?

Answer 1

Это то, что работает для меня. Возможно, можно использовать некоторую очистку.

class InterpolatedArray(object): 
    """ An array-like object that provides interpolated values between set points. 
    """ 
    points = None 
    wrap_value = None 
    offset = None 

    def _mod_delta(self, a, b): 
     """ Perform a difference within a modular domain. 
      Return a value in the range +/- wrap_value/2. 
     """ 
     limit = self.wrap_value/2. 
     val = a - b 
     if val < -limit: val += self.wrap_value 
     elif val > limit: val -= self.wrap_value 
     return val 

    def __init__(self, points, wrap_value=None): 
     """Initialization of InterpolatedArray instance. 

     Parameter 'points' is a list of two-element tuples, each of which maps 
     an input value to an output value. The list does not need to be sorted. 

     Optional parameter 'wrap_value' is used when the domain is closed, to 
     indicate that both the input and output domains wrap. For example, a 
     table of degree values would provide a 'wrap_value' of 360.0. 

     After sorting, a wrapped domain's output values must all be monotonic 
     in either the positive or negative direction. 

     For tables that don't wrap, attempts to interpolate values outside the 
     input range cause a ValueError exception. 
     """ 
     if wrap_value is None: 
      points.sort() # Sort in-place on first element of each tuple 
     else: # Force values to positive modular range 
      points = sorted([(p[0]%wrap_value, p[1]%wrap_value) for p in points]) 
      # Wrapped domains must be monotonic, positive or negative 
      monotonic = [points[x][1] < points[x+1][1] for x in xrange(0,len(points)-1)] 
      num_pos_steps = monotonic.count(True) 
      num_neg_steps = monotonic.count(False) 
      if num_pos_steps > 1 and num_neg_steps > 1: # Allow 1 wrap point 
       raise ValueError("Table for wrapped domains must be monotonic.") 
     self.wrap_value = wrap_value 
     # Pre-compute inter-value slopes 
     self.x_list, self.y_list = zip(*points) 
     if wrap_value is None: 
      intervals = zip(self.x_list, self.x_list[1:], self.y_list, self.y_list[1:]) 
      self.slopes = [(y2 - y1)/(x2 - x1) for x1, x2, y1, y2 in intervals] 
     else: # Create modular slopes, including wrap element 
      x_rot = list(self.x_list[1:]); x_rot.append(self.x_list[0]) 
      y_rot = list(self.y_list[1:]); y_rot.append(self.y_list[0]) 
      intervals = zip(self.x_list, x_rot, self.y_list, y_rot) 
      self.slopes = [self._mod_delta(y2, y1)/self._mod_delta(x2, x1) for x1, x2, y1, y2 in intervals] 

    def __getitem__(self, x):  # Works with indexing operator [] 
     result = None 
     if self.wrap_value is None: 
      if x < self.x_list[0] or x > self.x_list[-1]: 
       raise ValueError('Input value out-of-range: %s'%str(x)) 
      i = bisect.bisect_left(self.x_list, x) - 1 
      result = self.y_list[i] + self.slopes[i] * (x - self.x_list[i]) 
     else: 
      x %= self.wrap_value 
      i = bisect.bisect_left(self.x_list, x) - 1 
      result = self.y_list[i] + self.slopes[i] * self._mod_delta(x, self.x_list[i]) 
      result %= self.wrap_value 
     return result 

И тест:

import nose 

def xfrange(start, stop, step=1.): 
    """ Floating point equivalent to xrange().""" 
    while start < stop: 
     yield start 
     start += step 

# Test simple inverted mapping for non-wrapped domain 
pts = [(x,-x) for x in xfrange(1.,16., 1.)] 
a = InterpolatedArray(pts) 
for i in xfrange(1., 15., 0.1): 
    nose.tools.assert_almost_equal(a[i], -i) 
# Cause expected over/under range errors 
result = False # Assume failure 
try: x = a[0.5] 
except ValueError: result = True 
assert result 
result = False 
try: x = a[15.5] 
except ValueError: result = True 
assert result 

# Test simple wrapped domain 
wrap = 360. 
offset = 1.234 
pts = [(x,((wrap/2.) - x)) for x in xfrange(offset, wrap+offset, 10.)] 
a = InterpolatedArray(pts, wrap) 
for i in xfrange(0.5, wrap, 0.1): 
    nose.tools.assert_almost_equal(a[i], (((wrap/2.) - i)%wrap)) 

Answer 2

Но вы можете использование линейной интерполяции. Если ваше значение образца А равно, например, 7,75, что напоминает 2,5 градуса. Если значение образца B равно 180,35, оно также напоминает 2,5 градуса. Сложная часть - это когда переполнение значений, если это возможно вообще. Просто настройте кучу unittests, чтобы проверить, работает ли ваш алгоритм, и вы должны быстро идти.

Answer 3

Ответ на часть 1: таблица перевода, содержащая калибровочные значения + значение дрейфа.

В общем случае, если DialA сообщает 5.7, если физически оно находится на 0, 9.7, когда оно равно 5, тогда я установил бы значение дрейфа равным +/- .25 расстояния между каждой позицией считывания для учета механических и считывание дрейфа.

Ответ на часть 2: сохранение одинаковых значений на обоих циферблатах при отображении ожидаемой позиции.

Если вы не зависимы от направления, то просто вращайте выходной диск до тех пор, пока он не окажется в правильном положении в соответствии с вашей калибровочной таблицей.

Если вы зависимы от направления, вам нужно будет отслеживать последние 1-2 значения для определения направления. После того, как вы определили направление, вы можете переместить зависимый диск в нужном вам направлении, пока не будет достигнута позиция назначения.

Ваша калибровочная таблица также должна включать направление (положительное или отрицательное, например).

С учетом вышеуказанных двух частей вы сможете компенсировать смещения вращения и направленные переходы и производить точное считывание положения и направления.

Вот код, который дал таблицу калибровки, приведет к получению положения и направления, в котором будет решить проблему дисплея и делает матч зависимого набора с первичным набором:

#!/usr/bin/env python 

# Calibration table 
# calibrations[ device ][physical position]=recorded position 
calibrations={} 

calibrationsdrift=1.025 

calibrations["WheelA"]={} 

calibrations["WheelA"]={} 
calibrations["WheelA"]["000.0"]=5.7 
calibrations["WheelA"]["005.0"]=9.8 
calibrations["WheelA"]["010.0"]=13.9 
calibrations["WheelA"]["015.0"]=18.0 

calibrations["WheelB"]={} 
calibrations["WheelB"]["000.0"]=182.3 
calibrations["WheelB"]["005.0"]=178.4 
calibrations["WheelB"]["010.0"]=174.4 
calibrations["WheelB"]["015.0"]=170.3 


def physicalPosition(readout , device): 
     calibration=calibrations[device] 
     for physicalpos,calibratedpos in calibration.items(): 
       if readout < (calibratedpos + calibrationsdrift): 
         if readout > (calibratedpos - calibrationsdrift): 
           return physicalpos 
     return -0 


print physicalPosition(5.8 , "WheelA") 
print physicalPosition(9.8 , "WheelA") 
print physicalPosition(10.8 , "WheelA") 


def physicalDirection(lastreadout, currentreadout, device): 
     lastposition=physicalPosition(lastreadout, device) 
     currentposition=physicalPosition(currentreadout, device) 
     # Assumes 360 = 0, so 355 is the last position before 0 
     if lastposition < currentposition: 
       if lastposition == 000.0: 
         if currentposition == 355: 
           return -1 
       return 1 
     else: 
       return -1 


print physicalDirection(5.8, 10.8, "WheelA") 
print physicalDirection(10.8, 2.8, "WheelA") 

print physicalDirection(182, 174, "WheelB") 
print physicalDirection(170, 180, "WheelB") 

Запуском программы показывает что направление определяется правильно, даже для WheelB, который установлен назад на панель/устройства/и т.д.:

$ ./dials 
000.0 
005.0 
005.0 
1 
-1 
1 
-1 

Обратите внимание, что некоторые из значений «считывания», подаваемых в функции выключены. Это компенсируется значением дрейфа. Независимо от того, нужна ли вам одна, зависит от оборудования, с которым вы взаимодействуете.

Answer 4

Самое легкое решение состоит в том, чтобы все ваши элементы таблицы увеличивались (или уменьшались в зависимости от ситуации), добавляя или вычитая 360 для отдельных элементов, чтобы сделать это так. Двойной стол назад, чтобы он охватывал весь диапазон от 0 до 360 даже после всех дополнений и вычитаний. Это упрощает линейную интерполяцию. Затем вы можете взять по модулю 360 после расчета, чтобы вернуть его в диапазон.