Есть ли в наборах данных Era-interim неправильные строки proj4strings?

Я работаю с наборами данных Era-interim. Я хочу извлечь данные о погоде для некоторых городов. Код и данные обновляются до github.

Во-первых, я использую растр для чтения файла, который я скачал с сайта:

library(raster)
windspeed <- raster("data/10m_wind_speed_19950101.grib")
windspeed
# class       : RasterLayer 
# dimensions  : 241, 480, 115680  (nrow, ncol, ncell)
# resolution  : 0.75, 0.75  (x, y)
# extent      : -0.375, 359.625, -90.375, 90.375  (xmin, xmax, ymin, ymax)
# coord. ref. : +proj=longlat +a=6367470 +b=6367470 +no_defs 

Затем загружаю свои города:

load("capitals.RData")
head(capitals)
# ID iso3   country  capital   long    lat
# 1  1  AUS Australia Canberra 149.13 -35.31
# 2  2  AUT   Austria   Vienna  16.37  48.22
# 3  3  BEL   Belgium Brussels   4.33  50.83
# 4  4  BGR  Bulgaria    Sofia  23.31  42.69
# 5  5  BRA    Brazil Brasilia -47.91 -15.78
# 6  6  CAN    Canada   Ottawa -75.71  45.42

... и преобразовать их в объект sf:

library(sf)
capitals_sf <- st_as_sf(capitals, coords = c("long", "lat"), crs = 4326)
capitals_sf
# Simple feature collection with 40 features and 4 fields
# geometry type:  POINT
# dimension:      XY
# bbox:           xmin: -99.14 ymin: -35.31 xmax: 149.13 ymax: 60.17
# epsg (SRID):    4326
# proj4string:    +proj=longlat +datum=WGS84 +no_defs
# First 10 features:
#   ID iso3        country  capital              geometry
# 1   1  AUS      Australia Canberra POINT (149.13 -35.31)
# 2   2  AUT        Austria   Vienna   POINT (16.37 48.22)
# 3   3  BEL        Belgium Brussels    POINT (4.33 50.83)
# 4   4  BGR       Bulgaria    Sofia   POINT (23.31 42.69)
# 5   5  BRA         Brazil Brasilia POINT (-47.91 -15.78)
# 6   6  CAN         Canada   Ottawa  POINT (-75.71 45.42)
# 7   7  CHN          China  Beijing   POINT (116.4 39.93)
# 9   8  CYP         Cyprus  Nicosia   POINT (33.38 35.16)
# 11  9  CZE Czech Republic   Prague   POINT (14.43 50.08)
# 12 10  DEU        Germany   Berlin   POINT (13.38 52.52)

Поскольку windspeed и capital_sf имеют разные CRS, мне нужно выполнить некоторые преобразования:

newcrs <- crs(windspeed, asText=TRUE)
capitals_tf <- st_transform(capitals_sf, newcrs)
capital_tf
# Simple feature collection with 40 features and 4 fields
# geometry type:  POINT
# dimension:      XY
# bbox:           xmin: -99.14 ymin: -35.31 xmax: 149.13 ymax: 60.17
# epsg (SRID):    NA
# proj4string:    +proj=longlat +a=6367470 +b=6367470 +no_defs
# First 10 features:
#   ID iso3        country  capital              geometry
# 1   1  AUS      Australia Canberra POINT (149.13 -35.31)
# 2   2  AUT        Austria   Vienna   POINT (16.37 48.22)
# 3   3  BEL        Belgium Brussels    POINT (4.33 50.83)
# 4   4  BGR       Bulgaria    Sofia   POINT (23.31 42.69)
# 5   5  BRA         Brazil Brasilia POINT (-47.91 -15.78)
# 6   6  CAN         Canada   Ottawa  POINT (-75.71 45.42)
# 7   7  CHN          China  Beijing   POINT (116.4 39.93)
# 9   8  CYP         Cyprus  Nicosia   POINT (33.38 35.16)
# 11  9  CZE Czech Republic   Prague   POINT (14.43 50.08)
# 12 10  DEU        Germany   Berlin   POINT (13.38 52.52)

Странно, но proj4string меняется, но координаты не меняются.

Чтобы увидеть, удачно ли мое преобразование, я составляю сюжет:

plot(windspeed)
plot(capitals_tf, col = "black", add = TRUE)

вот сюжет:

Есть ли в наборах данных Era-interim неправильные строки proj4strings?

Долгота находится в диапазоне от -0,375 до 359,627, а не от -180 до 180. Таким образом, все города в восточном полушарии отмечены правильно, но все города в западном полушарии отсутствуют.

Я смущен. Почему не работает st_transform? Могу ли я передать неправильную строку proj4string, или функция просто не может обработать настроенную CRS?

Трудно помочь без воспроизводимый пример, но пока вы можете попробовать sp :: spTransform в качестве альтернативы, чтобы точно определить вашу проблему.

Michael Bird 10.08.2018 14:39

Спасибо большое. Пробовал sp :: spTransform. Это дало те же результаты. Я загрузил код и данные в github.com/yutuotuo84/era-interim-example. Он должен работать как воспроизводимый пример.

fishinweb 10.08.2018 17:25
Стоит ли изучать PHP в 2023-2024 годах?
Стоит ли изучать PHP в 2023-2024 годах?
Привет всем, сегодня я хочу высказать свои соображения по поводу вопроса, который я уже много раз получал в своем сообществе: "Стоит ли изучать PHP в...
Поведение ключевого слова "this" в стрелочной функции в сравнении с нормальной функцией
Поведение ключевого слова "this" в стрелочной функции в сравнении с нормальной функцией
В JavaScript одним из самых запутанных понятий является поведение ключевого слова "this" в стрелочной и обычной функциях.
Приемы CSS-макетирования - floats и Flexbox
Приемы CSS-макетирования - floats и Flexbox
Здравствуйте, друзья-студенты! Готовы совершенствовать свои навыки веб-дизайна? Сегодня в нашем путешествии мы рассмотрим приемы CSS-верстки - в...
Тестирование функциональных ngrx-эффектов в Angular 16 с помощью Jest
В системе управления состояниями ngrx, совместимой с Angular 16, появились функциональные эффекты. Это здорово и делает код определенно легче для...
Концепция локализации и ее применение в приложениях React ⚡️
Концепция локализации и ее применение в приложениях React ⚡️
Локализация - это процесс адаптации приложения к различным языкам и культурным требованиям. Это позволяет пользователям получить опыт, соответствующий...
Пользовательский скаляр GraphQL
Пользовательский скаляр GraphQL
Листовые узлы системы типов GraphQL называются скалярами. Достигнув скалярного типа, невозможно спуститься дальше по иерархии типов. Скалярный тип...
2
2
166
1

Ответы 1

Это хорошая справочная информация о формате набора данных ERA-Interim:

https://confluence.ecmwf.int/display/CKB/ERA-Interim%3A+What+is+the+spatial+reference

В нем конкретно говорится, что:

Longitudes range from 0 to 360, which is equivalent to -180 to +180 in Geographic coordinate systems.

Быстрый и грязный способ получить то, что вы хотите, - это "переместить" правую часть вашего растра. влево, затем вручную отрегулируйте экстент так, чтобы он составлял от -180 до 180. Таким образом, растр находится в "стандартном представлении GCS", и вы можете легко использовать его для черчения.

Например:

# create temporary raster, then "move" the right side to the left
tmp <- windspeed
tmp[, 1:240] <- windspeed[, 241:480]
tmp[, 241:480] <- windspeed[, 1:240]

# put data back in windspeed (not really needed) and update extent
windspeed <- tmp
extent(windspeed)@xmin <- extent(windspeed)@xmin -180
extent(windspeed)@xmax <- extent(windspeed)@xmax -180

windspeed
class       : RasterLayer 
dimensions  : 241, 480, 115680  (nrow, ncol, ncell)
resolution  : 0.75, 0.75  (x, y)
extent      : -180.375, 179.625, -90.375, 90.375  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +a=6367470 +b=6367470 +no_defs 
data source : in memory
names       : X10m_wind_speed_19950101 
values      : 0.9062432, 14.906  (min, max)

# now plot: 
capitals_sf <- st_as_sf(capitals, coords = c("long", "lat"), crs = 4326)

plot(windspeed)
plot(capitals_sf, col = "black", add = TRUE)

, что кажется более-менее правильным.

HTH!

Большое спасибо! Ваш ответ меня вдохновил!

fishinweb 11.08.2018 05:56

На самом деле я хочу обрабатывать множество файлов GRIB, используя один и тот же набор данных с заглавными буквами. Поэтому вместо того, чтобы менять местами растры, я поменял местами заглавные буквы: capitals <- capitals%>% mutate (long = ifelse (long <-0.375, long + 360, long)). Карта выглядит не очень хорошо, но функция извлечения дает тот же результат, что и я меняю растр местами.

fishinweb 11.08.2018 06:02

Другие вопросы по теме