КПД по данным телеметрии в Python

Интересный датасет нашел. Из области телеметрии. Набор данных (Data Set) подготовлен по результатам полетов устройства DJI Matrice 100. Ссылка в конце статьи. Данные представлены в виде временных рядов. Мы поработаем с отображением временных рядов и посчитаем энергозатраты и КПД по данным телеметрии.

Как данные попали в набор данных?

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

Во время каждого испытания изделия с уcтановленной ROS, генерируются файлы типа ROSbag. К набору данных разработчики приложили в том числе сырые файлы bag. Для нас сейчас не принципиально то, как именно набор bag файлов превращается в flights.csv. Мы просто отметим, что есть библиотеки для Python вроде Bagpy. Она ставится путем команды в активированном окружении:

pip install bagpy

Выведем ROSbag на печать (в консоль)

После импорта, можно импортировать и пользоваться данными из файла raw.bag, который мы для этого примера положили в корневую папку скрипта:

 import bagpy
from bagpy import bagreader
import pandas as pd

# Specify the path to your .bag file
bag_file = 'raw.bag'
# Create a bagreader object
bag = bagreader(bag_file)

# Print the list of topics available in the raw.csv file
print(&amp;amp;amp;amp;amp;amp;quot;Topics in the bag file:&amp;amp;amp;amp;amp;amp;quot;)
for topic in bag.topics:
    print(topic)

# Select the specific topic
anemometer_topic = '/anemometer'
# This will create csv file in raw subfolder called
# anemometer.csv and print full path.
csv_file = bag.message_by_topic(anemometer_topic)
print(csv_file)
# Read the CSV file into a DataFrame
df = pd.read_csv(csv_file)

# Print the first two messages
print(&amp;amp;amp;amp;amp;amp;quot;First message:&amp;amp;amp;amp;amp;amp;quot;)
print(df.iloc[0])
print(&amp;amp;amp;amp;amp;amp;quot;Second message:&amp;amp;amp;amp;amp;amp;quot;)
print(df.iloc[1])

Этот код выдаст в консоль все названия топиков, что есть в файле raw.bag и напечатает в консоль же первые два сообщения из топика ‘/anemometer’. Также он сохранит файл raw/anemometer.csv. Этот файл можно будет открывать блокнотом. Обратите внимание, что мы уже работаем в рамках Pandas стандартно загрузив в него данные из csv. Файл csv в этом скрипте был получен при помощи функции message_by_topic().

Построение графиков по при помощи Panel

Очень полезно строить графики, чтобы визуально оценивать величины. Мы будем строить графики при помощи Panel. Инструмент позволяет группировать графики в информационные панели и отображать их посредством доступа по сети к определенному порту через браузер.

Panel запускает сервер (адрес сообщает в консоли — например: http://localhost:5006), который показывает через браузер графики из Pyplot, Seaborn или Plotly. Если нужно порядок в раскладке, то графики желательно строить при помощи одного и того же инструмента. На рисунке выше видно, что интерактивный график выглядит маленьким.

В результате получаем следующие графики на одной веб-странице:

Из этого графика видно, что скорость имеет осцилляции в угловых точках траектории, где её вектор меняет направление. Скорость падает и снова набирает значение установки. Изделие останавливается (не до полной остановки) и набирает скорость в другом направлении.

По графику выше мы видим осцилляцию скорости в угловых точках траектории и устойчивую сходимость управления скоростью к значению установки (4 м/с). Время сходимости составляет около 15 секунд.

А это график траектории в проекции на плоскость земли. Север наверху. Красной точкой показано начало траектории. Изделие хорошо движется по прямой траектории. Чтобы не было каши из точек, мы сделали выборку каждой 10 точки из набора данных.

На графике траектории выше на осях не метры. Это широта и долгота. Поэтому север наверху. А шаги по значениям очень небольшие. Этот график в скрипте я поместил в комментарий, он не отображается при запуске скрипта.

Что же касается высоты, то мы можем наблюдать некоторый дрейф значения. Установленная высота полета равна 25 метрам. Однако, показания падают со временем. Это показано на рисунке ниже. Точность показаний датчика высоты — 0.5 метра.

Мы отображаем две оси одну слева — для Position Z, другую справа — для командного значения высоты. Видно, что командное значение всегда одинаковое, а фактическая высота падает.

Расчет энергозатрат

В рамках обработки телеметрии мы еще хотели бы оценить следующее. Нам примерно известна масса устройства от производителя:

2355 c батареей TB47D (4500 mAh применялась по Readme к data set)
2431 c батареей TB48D (5700 mAh применялась по статье в Nature)

На самом деле, в описании набора данных я указания о массе не нашел, однако есть модель изделия и есть ещё одна статья в журнале Nature (за авторством того же коллектива), где указано, что изделие функционирует с дополнительным оборудованием. Также, может быть дополнительная нагрузка. Поэтому, в файле flights.csv есть поле ‘payload’. Для хранения массы полетной нагрузки. Приведенные графики без нагрузки. Значит дополнительную взлетную массу дают датчики (анемометр) и их элементы крепления (анемометр устанавливается на штанге для невозмущенного потока).

Прибавим к массе производителя массу средств бортовых измерений (1.113 кг), получим: 2.431 + 1.113 + 0.136 = 3.680 (по статье в Nature).

Собственно КПД по данным телеметрии

Механика

Школьная физика, сначала Второй Закон Ньютона:

$F_{\text{total, 3D}} = a_{\text{total, 3D}} \times m$ где:

$F_{\text{total, 3D}}$ это равнодействующая сила (трехмерная, Н)

$a_{\text{total, 3D}}$ это полное ускорение (тоже трехмерное в м/с²)

Единицы измерения датасета даны в Readme. Там вообще много полезной информации.

Мощность механическую считаем так: $P_{\text{force}} = F_{\text{total, 3D}} \times V_{\text{total, 3D}}$ где:

$P_{\text{force}}$ это моментная мощность (Вт) по силе и скорости.

$V_{\text{total, 3D}}$ это моментная скорость (м/с).

Временные интервалы, на который мощность считается константой: $\Delta t = t_i - t_{i-1}$ .

Энергия на каждом временном шаге — это

$E_{\text{force}} = P_{\text{force}} \times \Delta t$

и полная произведенная механическая работа (энергия) это сумма:

$E_{\text{total, force}} = \sum E_{\text{force}}$

Электрическая энергия

С электрической стороны аналогично, вычисляем мощность, потом энергию по шагам и суммируя получаем полную:

$P_{\text{electric}} = V_{\text{battery}} \times I_{\text{battery}}$ где

$P_{\text{electric}}$ это электрическая мощность, далее — энергия на временных шагах:

$E_{\text{electric}} = P_{\text{electric}} \times \Delta t$

И сумма (то-есть общая затраченная энергия батареи):

$E_{\text{total, electric}} = \sum E_{\text{electric}}$

И теперь — интересно: эти два значения энергии можно сравнить и получить КПД. Вот результат работы приложенного ниже скрипта:

total_energy_force = 19440.235406235697 Дж
total_energy_electric = 78371.91157120022 Дж

По электроэнергии работа — почти 80 тысяч джоулей, а механическая менее 20 тысяч. Значит 25% КПД у модели. Это общий КПД (для этой конфигурации, траектории и погоды), в моментах времени он может сильно отличаться.

Я также построил график энергозатрат, но на таком отрезке времени это не информативно. Каша получается из пиков. Разве только видно, что в момент снижения, потребление по электрике то же, а механическая работа сильно меньше (КПД падает).

И ещё вот что очень важно здесь отметить. Мы в расчетах вообще не учли, что энергия ещё тратится на вращение, а угловое ускорение с моментами в расчете механической энергии не фигурировало. Это сделано по причине отсутствия данных об инерционных характеристиках конструкции. Вполне может быть, что вращение аппарата к существенным тратам энергии не относится. Я принял изделие за тонкое плоское тело с известной формулой моментов инерции. Посмотрел по телеметрии величины угловых ускорений пики там от 2 до 3 рад/с². В общем мощность на временном шаге получилась около 0.5 — 2.0 Вт. Я решил не вдаваться далее в расчеты, поскольку 19440 Дж / 150 сек = 126 Вт (средняя мощность линейной составляющей) — гораздо больше вращательной мощности.

Файлы для скачивания

Для работы скрипта в его папке должен быть файл flights.csv из датасета. Перечень необходимых библиотек определяется в секции импорта.

import pandas as pd
import seaborn as sns
import matplotlib.pyplot as plt
import panel as pn
import numpy as np
import plotly.graph_objects as go
from matplotlib.ticker import ScalarFormatter
import holoviews as hv

Скрипт отдельно:
Имя: UAV3DInteractiveTrajectory.7z
Размер: 2667 байтов (2 KiB)
CRC32: E4A89CCB
UAV3DInteractiveTrajectory Скачать

Датасет здесь.

Виртуальное окружение целиком — на Яндекс Диск.