Python-инструмент для получения данных о загрязнении воздуха из API качества воздуха Google Maps
Python-скрипт для получения данных о загрязнении воздуха через API Google Maps
Узнайте, как получить богатые данные о качестве воздуха в режиме реального времени со всего мира
В этой статье описано, как использовать API Google Maps Air Quality в Python для получения и исследования живых данных о загрязнении воздуха, временных рядов и карт. Проверьте полный кодздесь.
1. Основное
В августе 2023 года Google объявил о добавлении сервиса по качеству воздуха в список своих API для картографии. Вы можете прочитать больше об этом здесь. Кажется, эта информация теперь также доступна в приложении Google Maps, хотя данные, полученные с помощью API, оказались гораздо более полными.
Согласно объявлению, Google объединяет информацию из многих источников с различным разрешением – датчиков загрязнения на земле, спутниковых данных, информации о движении в режиме реального времени и прогнозов из численных моделей – для создания динамического набора данных о качестве воздуха в 100 странах с разрешением до 500 м. Звучит очень интересно и потенциально полезно для различных приложений в области картографии, здравоохранения и планирования!
Когда я впервые прочитал об этом, я планировал опробовать это в приложении “разговаривайте со своими данными”, используя некоторые вещи, изученные при создании этого инструмента для составления планов путешествий. Возможно, система, которая может строить временной ряд концентраций загрязнителей воздуха в вашем любимом городе, или, может быть, инструмент, помогающий людям планировать походы в своей местности, чтобы избежать плохого воздуха?
Здесь может помочь три инструмента API – сервис “текущие условия”, который предоставляет текущие значения индекса качества воздуха и концентрацию загрязняющих веществ в заданном месте; сервис “исторические условия”, который делает то же самое, но с интервалами в один час до 30 дней в прошлом, и сервис “тепловая карта”, который предоставляет текущие условия в заданной области в виде изображения.
Ранее я использовал отличный пакет googlemaps для вызова API Google Maps в Python, но эти новые API еще не поддерживаются. Удивительно, помимо официальной документации я нашел немного примеров использования этих новых инструментов и нет предварительно созданных пакетов Python, предназначенных для их вызова. Буду рад, если меня поправят, если кто-то знает другое!
Поэтому я создал свои собственные быстрые инструменты, и в этом посте мы рассмотрим, как они работают и как ими пользоваться. Надеюсь, это будет полезно для тех, кто хочет экспериментировать с этими новыми API в Python и ищет место для начала. Весь код для этого проекта можно найти здесь, и, вероятно, я буду расширять этот репозиторий со временем, добавляя больше функциональности и создавая некоторое приложение для картографии с данными о качестве воздуха.
2. Получение текущего качества воздуха в заданном месте
Давайте начнем! В этом разделе мы рассмотрим, как получить данные о качестве воздуха в заданном месте с помощью Google Maps. Вам сначала потребуется ключ API, который вы можете сгенерировать через свою учетную запись Google Cloud. У них есть бесплатный пробный период в течение 90 дней, после чего вы будете платить за использование API. Убедитесь, что вы включили “API качества воздуха” и ознакомьтесь с условиями ценообразования, прежде чем делать много запросов!
Снимок экрана библиотеки Google Cloud API, где вы можете активировать API для качества воздуха. Изображение, созданное автором.
Обычно я сохраняю свой ключ API в файле .env и загружаю его с помощью dotenv, используя функцию, подобную этой
Получение текущих условий требует POST-запроса, подробно описанного здесь. Мы возьмем вдохновение из пакета googlemaps, чтобы сделать это универсальным способом. Сначала мы создаем класс клиента, который использует requests для вызова. Цель довольно проста – мы хотим построить URL, аналогичный приведенному ниже, и включить все опции запроса, специфичные для запроса пользователя.
Класс Client принимает наш API-ключ как key и затем создает request_url для запроса. Он принимает параметры запроса в виде словаря params и помещает их в тело запроса JSON, обрабатываемое с помощью вызова self.session.post().
Теперь мы можем создать функцию, которая помогает пользователю собрать допустимые параметры запроса для API текущих условий, а затем использует этот класс Client для выполнения запроса. Снова, это вдохновлено дизайном пакета googlemaps.
def current_conditions(
client,
location,
include_local_AQI=True,
include_health_suggestion=False,
include_all_pollutants=True,
include_additional_pollutant_info=False,
include_dominent_pollutant_conc=True,
language=None,
):
"""
См. документацию для этого API здесь
https://developers.google.com/maps/documentation/air-quality/reference/rest/v1/currentConditions/lookup
"""
params = {}
if isinstance(location, dict):
params["location"] = location
else:
raise ValueError(
"Аргумент Location должен быть словарем, содержащим широту и долготу"
)
extra_computations = []
if include_local_AQI:
extra_computations.append("LOCAL_AQI")
if include_health_suggestion:
extra_computations.append("HEALTH_RECOMMENDATIONS")
if include_additional_pollutant_info:
extra_computations.append("POLLUTANT_ADDITIONAL_INFO")
if include_all_pollutants:
extra_computations.append("POLLUTANT_CONCENTRATION")
if include_dominent_pollutant_conc:
extra_computations.append("DOMINANT_POLLUTANT_CONCENTRATION")
if language:
params["language"] = language
params["extraComputations"] = extra_computations
return client.request_post("/v1/currentConditions:lookup", params)
Параметры для этого API довольно просты. Он требует словаря со значениями долготы и широты точки, которую вы хотите исследовать, и может опционально принимать различные другие аргументы, которые управляют объемом информации, возвращаемой в ответе. Давайте посмотрим его в действии с установленными всеми аргументами в значение True
Возвращается много интересной информации! У нас есть значения индекса качества воздуха от универсальных и американских индексов, а также концентрации основных загрязнителей, описание каждого из них и общий набор рекомендаций для здоровья в текущем состоянии воздуха.
{'dateTime': '2023-10-12T05:00:00Z', 'regionCode': 'us', 'indexes': [{'code': 'uaqi', 'displayName': 'Universal AQI', 'aqi': 60, 'aqiDisplay': '60', 'color': {'red': 0.75686276, 'green': 0.90588236, 'blue': 0.09803922}, 'category': 'Good air quality', 'dominantPollutant': 'pm10'}, {'code': 'usa_epa', 'displayName': 'AQI (US)', 'aqi': 39, 'aqiDisplay': '39', 'color': {'green': 0.89411765}, 'category': 'Good air quality', 'dominantPollutant': 'pm10'}], 'pollutants': [{'code': 'co', 'displayName': 'CO', 'fullName': 'Углекислый газ', 'concentration': {'value': 292.61, 'units': 'PARTS_PER_BILLION'}, 'additionalInfo': {'sources': 'Обычно возникает из-за неполного сгорания углеродсодержащих топлив, таких как в автомобильных двигателях и электростанциях.', 'effects': 'При вдыхании углекислого газа необходимое количество кислорода не поступает в кровь. Воздействие может вызывать головокружение, тошноту и головные боли. При высоких концентрациях возможна потеря сознания.'}}, {'code': 'no2', 'displayName': 'NO2', 'fullName': 'Азота диоксид', 'concentration': {'value': 22.3,
3. Получение временного ряда качества воздуха по заданному местоположению
Не было бы здорово иметь возможность получить временной ряд значений индекса качества воздуха (AQI) и загрязняющих веществ для заданного местоположения? Это могло бы раскрыть интересные закономерности, такие как корреляции между загрязняющими веществами или ежедневные флуктуации, вызванные движением транспорта или погодными факторами.
Мы можем сделать это при помощи другого POST-запроса к API "исторических условий", который даст нам почасовую историю. Это работает почти так же, как и текущие условия, единственное отличие заключается в том, что, поскольку результаты могут быть довольно длинными, они возвращаются в нескольких "страницах", что требует некоторой дополнительной логики для обработки.
Давайте изменить метод request_post класса Client для обработки этого.
def request_post(self,url,params): request_url = self.compose_url(url) request_header = self.compose_header() request_body = params response = self.session.post( request_url, headers=request_header, json=request_body, ) response_body = self.get_body(response) # поместим первую страницу в словарь ответа page = 1 final_response = { "page_{}".format(page) : response_body } # получим все страницы, если это необходимо while "nextPageToken" in response_body: # делаем вызов с токеном следующей страницы request_body.update({ "pageToken":response_body["nextPageToken"] }) response = self.session.post( request_url, headers=request_header, json=request_body, ) response_body = self.get_body(response) page += 1 final_response["page_{}".format(page)] = response_body return final_response
Этот код обрабатывает случай, когда в response_body содержится поле nextPageToken, которое является идентификатором следующей генерируемой и готовой к получению страницы данных. Если эта информация существует, нам просто нужно повторно вызвать API с новым параметром pageToken, который указывает на соответствующую страницу. Мы делаем это повторно в цикле while, пока не останется больше страниц. Наш словарь final_response теперь содержит еще один уровень, обозначенный номером страницы. Для вызовов current_conditions будет всегда только одна страница, но для вызовов historical_conditions их может быть несколько.
После выполнения этих задач мы можем написать функцию historical_conditions в очень похожем стиле, как для current_conditions.
def historical_conditions( client, location, specific_time=None, lag_time=None, specific_period=None, include_local_AQI=True, include_health_suggestion=False, include_all_pollutants=True, include_additional_pollutant_info=False, include_dominant_pollutant_conc=True, language=None,): """ См. документацию для данного API по ссылке https://developers.google.com/maps/documentation/air-quality/reference/rest/v1/history/lookup """ params = {} if isinstance(location, dict): params["location"] = location else: raise ValueError( "Аргумент location должен быть словарем, содержащим широту и долготу" ) if isinstance(specific_period, dict) and not specific_time and not lag_time: assert "startTime" in specific_period assert "endTime" in specific_period params["period"] = specific_period elif specific_time and not lag_time and not isinstance(specific_period, dict): # обратите внимание, что время должно быть в формате "Zulu" # например, datetime.datetime.strftime(datetime.datetime.now(),"%Y-%m-%dT%H:%M:%SZ") params["dateTime"] = specific_time # период задержки в часах elif lag_time and not specific_time and not isinstance(specific_period, dict): params["hours"] = lag_time else: raise ValueError( "Необходимо указать аргументы specific_time, specific_period или lag_time" ) extra_computations = [] if include_local_AQI: extra_computations.append("LOCAL_AQI") if include_health_suggestion: extra_computations.append("HEALTH_RECOMMENDATIONS") if include_additional_pollutant_info: extra_computations.append("POLLUTANT_ADDITIONAL_INFO") if include_all_pollutants: extra_computations.append("POLLUTANT_CONCENTRATION") if include_dominant_pollutant_conc: extra_computations.append("DOMINANT_POLLUTANT_CONCENTRATION") if language: params["language"] = language params["extraComputations"] = extra_computations # размер страницы по умолчанию установлен на 100 здесь params["pageSize"] = 100 # токен страницы будет заполнен, если это необходимо методом request_post params["pageToken"] = "" return client.request_post("/v1/history:lookup", params)
Для определения исторического периода API может принимать значение lag_time в часах, до 720 (30 дней). Он также может принимать словарь specific_period, который определяет начальное и конечное время в формате, описанном в комментариях выше. Наконец, для получения одного часа данных можно использовать только одну метку времени, предоставленную параметром specific_time. Также обратите внимание на использование параметра pageSize, который контролирует количество временных точек, возвращаемых в каждом вызове API. Значение по умолчанию здесь - 100.
Мы должны получить длинный вложенный JSON-ответ, который содержит значения индекса AQI и конкретные значения загрязняющих веществ с интервалом в 1 час за последние 720 часов. Существует много способов форматирования этого в структуру, удобную для визуализации и анализа, и функция ниже преобразует его в pandas dataframe в "long" формате, что хорошо работает с seaborn для построения графиков.
from itertools import chainimport pandas as pddef historical_conditions_to_df(response_dict): chained_pages = list(chain(*[response_dict[p]["hoursInfo"] for p in [*response_dict]])) all_indexes = [] all_pollutants = [] for i in range(len(chained_pages)): # необходимо проверить, если одна из меток времени не содержит данных, что иногда может случаться if "indexes" in chained_pages[i]: this_element = chained_pages[i] # получить время time = this_element["dateTime"] # получить все значения индекса и добавить метаданные all_indexes += [(time , x["code"],x["displayName"],"индекс",x["aqi"],None) for x in this_element['indexes']] # получить все значения загрязняющих веществ и добавить метаданные all_pollutants += [(time , x["code"],x["fullName"],"загрязняющее вещество",x["concentration"]["value"],x["concentration"]["units"]) for x in this_element['pollutants']] all_results = all_indexes + all_pollutants # сгенерировать dataframe в "long" формате res = pd.DataFrame(all_results,columns=["время","код","название","тип","значение","единицы"]) res["время"]=pd.to_datetime(res["время"]) return res
Запуск этого на выводе historical_conditions произведет dataframe, столбцы которого отформатированы для удобного анализа.
Универсальный AQI, AQI США, pm25 и pm10 значения для этого местоположения в ЛА на протяжении 30 дней. Изображение создано автором.
Это уже очень интересно! Очевидно, что во временных рядах загрязняющих веществ присутствуют несколько периодичностей, и важно отметить, что AQI США тесно коррелирует с концентрациями pm25 и pm10, как и ожидалось. Я менее знаком с Универсальным AQI, который предоставляет Google здесь, поэтому не могу объяснить, почему он антикоррелирован с pm25 и p10. Означает ли более низкий UAQI лучшее качество воздуха? Несмотря на некоторое поисковое исследование, мне так и не удалось найти хороший ответ.
4. Получение тепловых карт качества воздуха
Теперь перейдем к последнему случаю использования API Google Maps Air Quality - созданию тепловых карт. Документация по этой теме немного разрежена, что печально, потому что эти карты являются мощным инструментом для визуализации текущего качества воздуха, особенно в сочетании с картой Folium.
Мы получаем их с помощью GET-запроса, который включает построение URL в следующем формате, где местоположение тайла указывается с помощью zoom, x и y
GET https://airquality.googleapis.com/v1/mapTypes/{mapType}/heatmapTiles/{zoom}/{x}/{y}
Что означают zoom, x и y? Мы можем ответить на это, изучив, как Google Maps преобразует координаты широты и долготы в "координаты тайлов", что подробно описано здесь. По сути, Google Maps хранит изображения в сетках, где каждая ячейка имеет размеры 256 x 256 пикселей, и реальные размеры ячейки зависят от уровня масштабирования. При вызове API нам нужно указать, из какой сетки рисовать - это определяется уровнем масштабирования - и где на сетке рисовать - это определяется координатами тайлов x и y. В ответ получаем массив байтов, который может быть прочитан библиотекой обработки изображений на Python (PIL) или аналогичным пакетом.
После формирования нашего url по указанному формату, мы можем добавить несколько методов в класс Client, которые позволят нам получить соответствующее изображение.
def request_get(self,url): request_url = self.compose_url(url) response = self.session.get(request_url) # for images coming from the heatmap tiles service return self.get_image(response) @staticmethod def get_image(response): if response.status_code == 200: image_content = response.content # note use of Image from PIL here # needs from PIL import Image image = Image.open(io.BytesIO(image_content)) return image else: print("GET request for image returned an error") return None
Это хорошо, но нам действительно нужна возможность преобразовать набор координат широты и долготы в координаты тайлов. В документации объясняется, как это сделать - мы сначала преобразуем координаты в проекцию Меркатора, а затем преобразуем их в "пиксельные координаты", используя указанный уровень масштабирования. Наконец, мы преобразуем их в координаты тайлов. Для выполнения всех этих преобразований мы можем использовать класс TileHelper ниже.
import mathimport numpy as npclass TileHelper(object): def __init__(self, tile_size=256): self.tile_size = tile_size def location_to_tile_xy(self,location,zoom_level=4): # Основано на функции здесь # https://developers.google.com/maps/documentation/javascript/examples/map-coordinates#maps_map_coordinates-javascript lat = location["latitude"] lon = location["longitude"] world_coordinate = self._project(lat,lon) scale = 1 << zoom_level pixel_coord = (math.floor(world_coordinate[0]*scale), math.floor(world_coordinate[1]*scale)) tile_coord = (math.floor(world_coordinate[0]*scale/self.tile_size),math.floor(world_coordinate[1]*scale/self.tile_size)) return world_coordinate, pixel_coord, tile_coord def tile_to_bounding_box(self,tx,ty,zoom_level): # см. https://developers.google.com/maps/documentation/javascript/coordinates # для деталей box_north = self._tiletolat(ty,zoom_level) # номера тайлов развиваются к югу box_south = self._tiletolat(ty+1,zoom_level) box_west = self._tiletolon(tx,zoom_level) # номера тайлов развиваются к востоку box_east = self._tiletolon(tx+1,zoom_level) # (latmin, latmax, lonmin, lonmax) return (box_south, box_north, box_west, box_east) @staticmethod def _tiletolon(x,zoom): return x / math.pow(2.0,zoom) * 360.0 - 180.0 @staticmethod def _tiletolat(y,zoom): n = math.pi - (2.0 * math.pi * y)/math.pow(2.0,zoom) return math.atan(math.sinh(n))*(180.0/math.pi) def _project(self,lat,lon): siny = math.sin(lat*math.pi/180.0) siny = min(max(siny,-0.9999), 0.9999) return (self.tile_size*(0.5 + lon/360), self.tile_size*(0.5 - math.log((1 + siny) / (1 - siny)) / (4 * math.pi))) @staticmethod def find_nearest_corner(location,bounds): corner_lat_idx = np.argmin([ np.abs(bounds[0]-location["latitude"]), np.abs(bounds[1]-location["latitude"]) ]) corner_lon_idx = np.argmin([ np.abs(bounds[2]-location["longitude"]), np.abs(bounds[3]-location["longitude"]) ]) if (corner_lat_idx == 0) and (corner_lon_idx == 0): # ближайший - latmin, lonmin direction = "юго-запад" elif (corner_lat_idx == 0) and (corner_lon_idx == 1): direction = "юго-восток" elif (corner_lat_idx == 1) and (corner_lon_idx == 0): direction = "северо-запад" else: direction = "северо-восток" corner_coords = (bounds[corner_lat_idx],bounds[corner_lon_idx+2]) return corner_coords, direction @staticmethod def get_ajoining_tiles(tx,ty,direction): if direction == "юго-запад": return [(tx-1,ty),(tx-1,ty+1),(tx,ty+1)] elif direction == "юго-восток": return [(tx+1,ty),(tx+1,ty-1),(tx,ty-1)] elif direction == "северо-запад": return [(tx-1,ty-1),(tx-1,ty),(tx,ty-1)] else: return [(tx+1,ty-1),(tx+1,ty),(tx,ty-1)]
Мы видим, что location_to_tile_xy получает словарь местоположения и уровень масштабирования, а затем возвращает тайл, в котором можно найти эту точку. Другая полезная функция - tile_to_bounding_box, которая найдет охватывающие координаты указанной сетки. Нам это понадобится, если мы собираемся определить местоположение ячейки и отобразить ее на карте.
Давайте посмотрим, как это работает внутри функции air_quality_tile ниже, которая будет принимать наш client, location и строку, указывающую, какой тип тайла мы хотим получить. Нам также нужно указать уровень масштабирования, выбор которого может быть сложным в начале и требует некоторой проб и ошибок. Мы обсудим аргумент get_adjoining_tiles в скором времени.
def air_quality_tile( client, location, pollutant="UAQI_INDIGO_PERSIAN", zoom=4, get_adjoining_tiles = True): # см. https://developers.google.com/maps/documentation/air-quality/reference/rest/v1/mapTypes.heatmapTiles/lookupHeatmapTile assert pollutant in [ "UAQI_INDIGO_PERSIAN", "UAQI_RED_GREEN", "PM25_INDIGO_PERSIAN", "GBR_DEFRA", "DEU_UBA", "CAN_EC", "FRA_ATMO", "US_AQI" ] # содержит полезные методы для работы с координатами тайла helper = TileHelper() # получаем тайл, в котором находится местоположение world_coordinate, pixel_coord, tile_coord = helper.location_to_tile_xy(location,zoom_level=zoom) # получаем ограничивающий прямоугольник тайла bounding_box = helper.tile_to_bounding_box(tx=tile_coord[0],ty=tile_coord[1],zoom_level=zoom) if get_adjoining_tiles: nearest_corner, nearest_corner_direction = helper.find_nearest_corner(location, bounding_box) adjoining_tiles = helper.get_ajoining_tiles(tile_coord[0],tile_coord[1],nearest_corner_direction) else: adjoining_tiles = [] tiles = [] # получаем все смежные тайлы, плюс тот, который нас интересует for tile in adjoining_tiles + [tile_coord]: bounding_box = helper.tile_to_bounding_box(tx=tile[0],ty=tile[1],zoom_level=zoom) image_response = client.request_get( "/v1/mapTypes/" + pollutant + "/heatmapTiles/" + str(zoom) + '/' + str(tile[0]) + '/' + str(tile[1]) ) # преобразуем изображение PIL в массив numpy try: image_response = np.array(image_response) except: image_response = None tiles.append({ "bounds":bounding_box, "image":image_response }) return tiles
Из чтения кода мы видим, что рабочий процесс выглядит следующим образом: сначала находим координаты тайла интересующего нас местоположения. Это определяет сеточную ячейку, которую мы хотим получить. Затем находим ограничивающие координаты этой сеточной ячейки. Если мы хотим получить окружающие тайлы, находим ближайший угол ограничивающего прямоугольника и затем используем его для вычисления координат тайлов трех смежных сеточных ячеек. Затем вызываем API и возвращаем каждый из тайлов как изображение с соответствующим ограничивающим прямоугольником.
Мы можем запустить это стандартным образом, следующим образом:
А затем отобразить с помощью folium для карты с возможностью масштабирования! Обратите внимание, что здесь я использую leafmap, потому что пакет может создавать карты Folium, совместимые с gradio, мощным инструментом для создания простых пользовательских интерфейсов для приложений на Python. Взгляните на эту статью для примера.
Возможно, разочаровывающим образом, на этом уровне масштабирования плитка, содержащая наше местоположение, в основном состоит из моря, хотя все же интересно видеть отображение загрязнения воздуха на подробной карте. Если вы увеличите масштаб, можно увидеть, что информация о дорожном движении используется для информирования сигналов о качестве воздуха в городских районах.
Отображение плитки тепловой карты качества воздуха на верхушке карты Folium. Изображение сгенерировано автором.
Установка get_adjoining_tiles=True дает нам намного более приятную карту, поскольку она извлекает три ближайшие не перекрывающиеся плитки на этом уровне масштабирования. В нашем случае это существенно помогает сделать карту более презентабельной.
Когда мы также извлекаем смежные плитки, получается более интересный результат. Заметьте, что здесь показаны универсальные показатели качества воздуха AQI. Изображение сгенерировано автором.
Я лично предпочитаю изображения, сгенерированные при pollutant=US_AQI, но есть несколько различных вариантов. К сожалению, API не возвращает цветовую шкалу, хотя было бы возможно создать ее, используя значения пикселей в изображении и знание о том, что означают цвета.
Те же плитки, что и выше, окрашены согласно стандарту AQI США. Эта карта была сгенерирована 10/12/2023 года, и яркое пятно красного цвета в центральной части штата Калифорния, по-видимому, вызвано контролируемым сжиганием на холмах рядом с Коалингой, согласно этому инструменту https://www.frontlinewildfire.com/california-wildfire-map/. Изображение сгенерировано автором.
Заключение
Спасибо, что дочитали до конца! Здесь мы исследовали, как использовать API качества воздуха Google Maps для получения результатов в Python, которые могут быть использованы в интересных приложениях. В будущем я надеюсь продолжить другой статьей о инструменте air_quality_mapper, поскольку он дальше развивается, но надеюсь, что здесь обсуждаемые скрипты также окажутся полезными по своему назначению. Как всегда, буду рад любым предложениям по дальнейшему улучшению!