[Перевод] Анализ данных погодной станции, основанной на Arduino13.05.2015 13:17

Создание собственной личной метеостанции стало намного проще, чем раньше. С учетом непостоянной погоды в Новой Англии, мы решили, что хотим создать нашу собственную метеостанцию и использовать MATLAB для анализа метео данных.В статье мы ответим на следующие вопросы:

В каком направлении дул ветер в течение последних 3-х часов? Как изменялись температура и точка росы в течение последней недели? На самом ли деле падает барометрическое давление при приближении грозы? Понятно, что рассмотренные вопросы достаточно просты, но описанные приемы и команды помогут вам решать более сложные практические задачи.

Шаг 1: Размещение метеостанцииВо-первых, мы должны были решить, где разместить нашу метеостанцию. Мы решили сделать это на верхнем этаже автостоянки. Место было выбрано, поскольку оно было под воздействием погодных явлений, но также и имело крышу для монтажа электроники. Из-за того, что мы в конечном счете хотели отдавать данные стороннему агрегатору данных, мы должны были учитывать это при выборе места и железа.574408735637416daa6802671e48e7b6.pngШаг 2: Выбор аппаратного обеспечения Выбранное расположение было вне радиуса действия Wi-Fi нашего здания, поэтому мы должны были найти способ передачи данных от станции к приемнику, который находился внутри соседнего здания.Чтобы сделать это, мы подсоединили Arduino® Uno к передатчику высокой мощности XBee®. Затем данные передаются в принимающий модуль XBee на Arduino Mega уже внутри здания. Этот Arduino был подключен к Интернету и данные отправляются в бесплатный сервис агрегации данных, ThingSpeak.com. Полный список аппаратуры, которую мы использовали в нашей установке приведен ниже.Список оборудования 44162227778743888ac368fd68d19a5b.pngШаг 3: Подключение передатчика метеостанции и программирование наружной Arduino Внешнее Arduino установлено на гараже и отвечает за сбор измерений и передачу данных на Arduino в помещении. Для реализации такого обмена, сначала мы соединили контакты от Weather shield и XBee shield. Затем мы соединили XBee shield с Arduino Uno в соответствии с документацией, предоставленной с шилдами. В XBee shield есть XBee трансивер высокой мощности. Мы использовали программное обеспечение X-CTU для программирования нужного адреса получателя и загрузки нужной прошивки в XBee shield. Анемометр, датчики ветра дождя были подключены к weather shield с помощью разъемов RJ-45, входящими в комплект поставки.

Шаг 4: Подключение приемника погодной станции и программирование Arduino, находящегося в помещении Внутреннее Arduino находится внутри нашего здания, и отвечает за получение данных от наружного Arduino, за проверку достоверности данных, и их отправку на ThingSpeak.Чтобы построить такую систему, сначала мы припаяли контакы на Ethernet shield и второй XBee shield.Затем мы вставили Ethernet shield и XBee shield в Мега Arduino в соответствии с документацией. Мы вновь использовали X-CTU для программирования необходимого адреса получателя и загрузки нужной прошивку на XBee shield. Далее, мы запрограммировали Arduino получить XBee сообщения и пересылать пакеты на ThingSpeak со скоростью один раз в минуту. Перед отправкой сообщения ThingSpeak, мы настроили учетную запись и сконфигурировали информацию о канале и расположении. Как только закончилась настройка передачи данных, мы начали анализировать данные в MATLAB.Обратите внимание, что для воспроизведения анализа, приведенного в статье, вам не нужно покупать и настраивать оборудование.Текущие данные нашей станции доступны на канале 12 397.

Шаг 5: Отвечаем на вопросы Получение данных о погоде из ThingSpeak Чтобы ответить на наши первые два вопроса, мы используем команду thingSpeakFetch для просмотра доступных полей данных, одновременного импорта всех полей в MATLAB и хранения данных о температуре, влажности, скорости ветра, и направления ветра в своих переменных. Документация о поддержке ThingSpeak. [d, t, ci] = thingSpeakFetch (12397,'NumPoints',8000); % fetch last 8000 minutes of data 8000 точек — это максимальное количество точек, которое ThingSpeak позволяет запросить за раз. Для нашей частоты измерений, это соответствует примерно 6 дням измерений. tempF = d (:,4); % field 4 is temperature in deg F baro = d (:,6); % pressure in inches Hg humidity = d (:,3); % field 3 is relative humidity in percent windDir = d (:,1); windSpeed = d (:,2); tempC = (5/9)*(tempF-32); % convert to Celsius availableFields = ci.FieldDescriptions' availableFields =

'Wind Direction (North = 0 degrees)' 'Wind Speed (mph)' '% Humidity' 'Temperature (F)' 'Rain (Inches/minute)' 'Pressure («Hg)' 'Power Level (V)' 'Light Intensity' Вычисление некоторых основных статистических данных за период наблюдений и их визуализация Чтобы получить лучшее понимание наших данных, мы находим минимальное, максимальное и среднее значения для данных, которые мы импортировали и мы находим время, для максимальных и минимальных значений.Это дает нам быстрый путь к валидации данных с нашей метеостанции. [maxData, index_max] = max (d); maxData = maxData'; times_max = datestr (t (index_max)); times_max = cellstr (times_max); [minData, index_min] = min (d); minData = minData'; times_min = datestr (t (index_min)); times_min = cellstr (times_min); meanData = mean (d); meanData = meanData'; % make column vector summary = table (availableFields, maxData, times_max, meanData, minData, times_min) % display summary =

availableFields maxData times_max ____________________________________ _______ ______________________

'Wind Direction (North = 0 degrees)' 338 '10-Jul-2014 05:01:32' 'Wind Speed (mph)' 6.3 '10-Jul-2014 12:47:14' '% Humidity' 86.5 '15-Jul-2014 04:51:24' 'Temperature (F)' 96.7 '12-Jul-2014 16:28:55' 'Rain (Inches/minute)' 0.04 '15-Jul-2014 13:47:13' 'Pressure («Hg)' 30.23 '11-Jul-2014 09:25:07' 'Power Level (V)' 4.44 '10-Jul-2014 10:25:01' 'Light Intensity' 0.06 '12-Jul-2014 13:23:38'

meanData minData times_min _________ _______ ______________________

NaN 0 '10-Jul-2014 04:54:32' 3.0272 0 '10-Jul-2014 01:33:14' 57.386 25.9 '12-Jul-2014 13:39:39' 80.339 69.6 '11-Jul-2014 06:59:54' 5.625e-05 0 '10-Jul-2014 01:02:11' 30.04 29.78 '15-Jul-2014 13:04:08' 4.4149 4.38 '11-Jul-2014 09:22:06' 0.0092475 0 '10-Jul-2014 01:02:11' Если мы получаем неожиданные значения, такие как максимальное атмосферное давление 40 атмосфер или максимальную температуру 1700 градусов, мы можем предположить, что данные некорректны. Такие ошибки могут происходить из-за ошибок передачи, выбросов напряжения и других причин. В конце статьи мы покажем некоторые способы обработки таких «выбросов», но для выгруженных данных, когда этот отчет был опубликован, все выглядит в порядке.Для приведенного выше кода, мы использовали табличный тип данных, впервые введенный в MATLAB R2013b.См запись в блоге для дальнейшего ознакомления с этим типом данных.Визуализация ветра за последние три часа Так как наша погодная станция получает метеорологические сводки примерно раз в минуту, мы смотрим на последние 180 минут, что бы ответить на наш вопрос о том, какой ветер был в последние 3 часа, и мы используем Compass plot в MATLAB, чтобы увидеть скорость и направление ветра в интересующий времени. Это математический компас, где Север (0 градусов) находится справа, и значения градуса увеличиваются против часовой стрелки. 90 градусов представляет Восток (вверху), 180 градусов представляет юг (слева) и 270 градусов представляет Запад (внизу). Если не произошло сдвига ветра из-за грозы или при лобовом проходе, компас, как правило, показывает преимущественное направление ветра, обозначенное стрелками более высокой плотности на графике. В случае данных, запрошенных в этом примере мы наблюдаем направление ветра преимущественно с юга и юго-запада, типичное для летнего времени в Новой Англии. figure (1) windDir = windDir ((end-180): end); % last 3 hours windSpeed = windSpeed ((end-180): end); rad = windDir*2*pi/360; u = cos (rad) .* windSpeed; % x coordinate of wind speed on circular plot v = sin (rad) .* windSpeed; % y coordinate of wind speed on circular plot compass (u, v) titletext = strcat (' for last 3 hours ending on: ', datestr (now)); title ({'Wind Speed, MathWorks Weather Station'; titletext}) d1acb57c7dbf41fcb3a94288bf83664e.pngРасчет точки росы Теперь мы готовы ответить на наш второй вопрос о том, как температура и точка росы менялись в течениипрошлой недели.Точка росы — это температура, при которой воздух (при охлаждении) становится насыщенным водяным паром. Чем более влажность воздушной массы, тем выше точка росы. Точка росы также иногда используется как мера дискомфорта. Когда точка росы превышает 65 градусов (+18С), многие люди начинают говорить, что воздух «липкий». При точке росы за 70 многие чувствуют себя некомфортно. Общая оценка для точки росы, TDP, может быть найдена с помощью уравнения и постоянных, как показано ниже: f6c8f5f70d3a43fbadac9c177a5699bc.pngгде098d82c189ef4a8db105c104abb6da14.pngс 1199eaeaa98d412fb461df84a013dcb4.pngТеперь, написав приведенные выше уравнения, как код MATLAB мы имеем: b = 17.62; c = 243.5; gamma = log (humidity/100) + b*tempC ./ (c+tempC); tdp = c*gamma ./ (b-gamma); tdpf = (tdp*1.8) + 32; % convert back to Fahrenheit Визуализация температуры, влажности воздуха и точки росы в зависимости от времени Теперь, когда мы рассчитали точку росы, мы готовы визуализировать данные и наблюдать их поведение в течение последних 5 или 6 дней. figure (2) [ax, h1, h2] = plotyy (t,[tempF tdpf], t, humidity); set (ax (2),'XTick',[]) set (ax (2),'YColor','k') set (ax (1),'YTick',[0,20,40,60,80,100]) set (ax (2),'YTick',[0,20,40,60,80,100]) datetick (ax (1),'x','keeplimits','keepticks') set (get (ax (2),'YLabel'),'String', availableFields (3)) set (get (ax (1),'YLabel'),'String', availableFields (4)) grid on legend ('Location','South','Temperature','Dew point', 'Humidity') 7413587157124fdabd73d18601592434.pngВо время обычной недели, вы можете ясно видеть суточные изменения температуры и влажности. Как и ожидалось относительная влажность, как правило, возрастает ночью (при понижении температуры к точке росы) и максимальные температуры, как правило, достигаются во второй половине дня.Температура точки росы указывает влажность воздушных масс. Когда мы собрали данные для публикации этого примера, вы можете видеть, что точка росы была более 70 градусов, что является типичной для жаркого и влажного летнего дня в Новой Англии.Если вы запустите этот код в MATLAB, вы можете получить разные ответы, по мере того, как вы будете запрашивать самые последние данные, сообщенные метеостанцией.Получение и обработка данных о погоде за дождливый день в Новой Англии Последний вопрос, на который мы хотели ответить — действительно ли атмосферное давление падает перед дождем? Чтобы сделать это, мы получили данные из метеостанции за известный дождливый день. На этот раз, мы заинтересованы в барометрическом давлении и количестве осадков. Мы использовали датчик самостоятельного опорожнения, который опорожняется, когда он заполнен. Наш датчик вращается и опорожняется, когда выпало 0,01 дюймов осадков.Наш код Arduino считает количество опорожненийза каждую минуту и передает соответствующее значение осадков в ThingSpeak. Мы используем MATLAB для выборки из наших данных почасовых образцов, чтобы мы могли легко увидеть накопление осадков и тенденции давления. [d, t, ci] = thingSpeakFetch (12397,'DateRange',{'6/4/2014','6/6/2014'}); % get data baro = d (:,6); % pressure extraData = rem (length (baro),60); % computes excess points beyond the hour baro (1: extraData) = []; % removes excess points so we have even number of hours rain = d (:,5); % rainfall from sensor in inches per minute Действительно ли шел ливень 5 июня? Это трудно сказать, если мы просто посмотрим на каждую минуту передачи данных, так как максимальное значение составляет всего 0,01 дюйма. Однако, если суммировать все осадки за 5 июня, мы видим, что мы получили 0,48 дюйма осадков, что составляет 13% от среднемесячного 3,68 дюйма, что показывает, что день был действительно очень дождливым. Чтобы получить более полное представление о том, когда выпал максимум осадков, мы преобразовали данные в часовые данные, как показано ниже. rain (1: extraData) = []; t (1: extraData) = []; rainHourly = sum (reshape (rain,60,[]))'; % convert to hourly summed samples maxRainPerMinute = max (rain) june5rainfall = sum (rainHourly (25: end)) % 24 hours of measurements from June 5 baroHourly = downsample (baro,60); % hourly samples timestamps = downsample (t,60); % hourly samples maxRainPerMinute =

0.0100

june5rainfall =

0.4800 Визуализация данных о почасовой облачности и давления и нахождение тренда давления После того как мы предварительно обработали наши данные, мы готовы визуализировать их. Здесь мы используем MATLAB, для нахождения линии тренда к данным барометрического давления. Если построить график, увидим ли мы падение давления перед тяжелым ливнем? figure (3) subplot (2,1,1) bar (timestamps, rainHourly) % plot rain xlabel ('Date and Time') ylabel ('Rainfall (inches /per hour)') grid on datetick ('x','dd-mmm HH: MM','keeplimits','keepticks') title ('Rainfall on June 4 and 5') subplot (2,1,2) hold on plot (timestamps, baroHourly) % plot barometer xlabel ('Date and Time') ylabel (availableFields (6)) grid on datetick ('x','dd-mmm HH: MM','keeplimits','keepticks') detrended_Baro = detrend (baroHourly); baroTrend = baroHourly — detrended_Baro; plot (timestamps, baroTrend,'r') % plot trend hold off legend ('Barometric Pressure','Pressure Trend') title ('Barometric Pressure on June 4 and 5') 48a80aa664304ec0bd1a31b65b074e97.pngПосле того, как мы визуализируем эти данные и посмотрим на тренд, мы ясно видим, что атмосферное давление действительно падает ровно до начала осадков! Очистка данных от недостоверных измерений Наши вопросы были довольно простыми, но иногда возникают проблемы с имеющимися данными. Например, 30 мая, мы записали несколько недостоверных данных о температуре. Давайте посмотрим, как использовать MATLAB для их фильтрации. [d, t] = thingSpeakFetch (12397,'DateRange',{'5/30/2014','5/31/2014'}); rawTemperatureData = d (:,4); newTemperatureData = rawTemperatureData; minTemp = min (rawTemperatureData) % wow that is cold! minTemp =

-1.7662e+03 Использование порогового фильтра для удаления выброса данных Удаление элементов, которые не прошли пороговый тест. В этом случае, у нас есть некоторые значения, которые заведомо недостоверны, такие как значение температуры -1766 градусов по Фаренгейту. Поэтому можно использовать данные, которые включают в себя только значения температуры, между 0 и 120, которые являются приемлемыми величинами для весеннего дня в Новой Англии. tnew = t'; outlierIndices = [find (rawTemperatureData < 0); find(rawTemperatureData > 120)]; tnew (outlierIndices) = []; newTemperatureData (outlierIndices) = []; Построим графики очищенных и исходных данных. figure (4) subplot (2,1,2) plot (tnew, newTemperatureData,'-o') datetick xlabel ('Time of Day') ylabel (availableFields (4)) title ('Filtered Data — outliers deleted') grid on subplot (2,1,1) plot (t, rawTemperatureData,'-o') datetick xlabel ('Time of Day') ylabel (availableFields (4)) title ('Original Data') grid on 673bd527ee0c4f2e8429147f3b9a0d4b.pngИспользование медианного фильтра для удаления недостоверных данных Еще один способ удаления плохих данных — применить медианный фильтр. Медианный фильтр не требует столько знаний о наборе данных. Он будет просто удаляет значения, которые оказываются вне среднего ближайших соседей. Результаты фильтрации — вектор той же длины что и исходный, в отличие от удаление точек данных, что приводит к разрывам в данных и укорачиванию записи. Этот тип фильтра также может быть использован для удаления шума из сигнала. n = 5; % this value determines the number of total points used in the filter Большие значения n обозначают количество «соседей» для сравнения. При температурах, собранных раз в минуту, мы выбираем n = 5, потому что температура не должна, как правило, сильно меняться в течение 5 минут. f = medfilt1(rawTemperatureData, n); figure (5) subplot (2,1,2) plot (t, f,'-o') datetick xlabel ('Time of Day') ylabel (availableFields (4)) title ('Filtered Data — using median filter') grid on subplot (2,1,1) plot (t, rawTemperatureData,'-o') datetick xlabel ('Time of Day') ylabel (availableFields (4)) title ('Original Data') grid on cbd321a662184f60b85ea1af93d3bea3.pngВывод Итак, мы проанализировали данные с нашей метеостанции и ответили на интересующие нас вопросы.Мы также продемонстрировали некоторые способы фильтрации данных, которые не проходят первоначальную проверку.Вопросы, поставленные в статье, элементарны, но подумайте чего можно добиться собирая данные с тысячи станций? Кроме этого можно измерять и анализировать real-time изменения промышленных и научных систем и принимать решения на основе данных.

Напоминаем, что для воспроизведения анализа, приведенного в статье, вам не нужно покупать и настраивать оборудование. Оперативные данные из нашей установки доступны на канале 12 397 на ThingSpeak.com.

MATLAB код проекта

© Habrahabr.ru