Преобразование экстента растрового стека из метров в координаты широты / долготы
У меня есть растровый стек в следующем формате, который я преобразовал в NetCDF, используя этот метод. Это работает, но переменные широты и долготы находятся в "метрах", что является экстентом растрового файла, и мне нужно, чтобы они были в десятичных градусах. Файл R Grid находится здесь и имеет следующий формат:
class : RasterBrick
dimensions : 205, 170, 34850, 12 (nrow, ncol, ncell, nlayers)
resolution : 1, 1 (x, y)
extent : 160.5, 330.5, 145.5, 350.5 (xmin, xmax, ymin, ymax)
coord. ref. : NA
data source : \TavgM_1981.grd
names : Jan.1981, Feb.1981, Mar.1981, Apr.1981, May.1981, Jun.1981, Jul.1981, Aug.1981, Sep.1981, Oct.1981, Nov.1981, Dec.1981
min values : 1.9912137, 0.8120775, 4.0446638, 4.4135274, 6.5349769, 8.6149150, 9.9991916, 11.8400562, 9.6407796, 3.5005649, 4.1205872, -0.6106244
max values : 9.850221, 9.121176, 10.238524, 9.858942, 11.669445, 14.260988, 15.722292, 17.235090, 15.708690, 11.598482, 11.552235, 8.981533
time : Jan 1981, Feb 1981, Mar 1981, Apr 1981, May 1981, Jun 1981, Jul 1981, Aug 1981, Sep 1981, Oct 1981, Nov 1981, Dec 1981
В идеале я хотел бы преобразовать это в файл NetCDF, который имеет следующий формат для чтения моделью:
Dimensions: (lat: 408, lon: 881, nb2: 2, time: 660)
Coordinates:
* lon (lon) float64 -44.81 -44.69 -44.56 -44.44 -44.31 -44.19 ...
* lat (lat) float64 21.81 21.94 22.06 22.19 22.31 22.44 22.56 22.69 ...
* time (time) datetime64[ns] 1951-01-16T12:00:00 1951-02-16T12:00:00 ...
Dimensions without coordinates: nb2
Data variables:
time_bnds (time, nb2) datetime64[ns] 1951-01-16T12:00:00 ...
tas (time, lat, lon) float64 nan nan nan nan nan nan nan nan nan ...
Я попытался преобразовать экстент Растрового кирпича в широты и долготы, используя этот метод:
sputm <- SpatialPoints(test@extent@xmin, proj4string=CRS("+proj=utm +zone=29N +datum=WGS84"))
Но это дает эту ошибку:
Error in (function (classes, fdef, mtable) : unable to find an inherited method for function ‘coordinates’ for signature ‘"numeric"’
РЕДАКТИРОВАТЬ: также пытался
crs(test) <- '+proj=merc +datum=WGS84'
x <- projectRaster(test, crs='+proj=longlat +datum=WGS84')
writeRaster(x, "rstack2.nc", overwrite=TRUE, format="CDF", varname="Temperature", varunit="degC",
longname="Temperature -- raster stack to netCDF", xname="Longitude", yname="Latitude", zname="Time (Month)")
Но когда я рассматриваю это в общем, координаты широты и долготы равны 0,0.
Я не уверен, куда идти отсюда.
1 ответ
Как только я узнал Ирландию, я посмотрел на прогнозы Ирландии. Я думаю, что это TM75 / Irish Grid. Я изменил разрешение вашего растра так, чтобы оно составляло 2500м.
library(raster)
r <- stack("~/Bureau/TavgM_1981")
xmax(r) <- xmin(r) + 2500 * ncol(r)
ymax(r) <- ymin(r) + 2500 * nrow(r)
crs(r) <- "+proj=tmerc +lat_0=53.5 +lon_0=-8 +k=1.000035 +x_0=200000 +y_0=250000 +ellps=mod_airy +towgs84=482.5,-130.6,564.6,-1.042,-0.214,-0.631,8.15 +units=m +no_defs"
# Convert to degrees
x_wgs84 <- projectRaster(r, crs='+proj=longlat +datum=WGS84')
writeRaster(x, "~/Bureau/TavgM_wgs84.tif", overwrite = TRUE)
Ниже вы увидите результат этого разрешения + изменение проекции по сравнению с целью.
Мы не так далеки от реальности, но в вашем оригинальном растре есть что-то плохое. Если это (к сожалению) обычно не иметь crs, когда кто-то дает вам слой, это обычно не пропустить правильное разрешение.
Поскольку предоставленный вами растр сохраняется как файл ".gri /.grd", это означает, что он был создан с использованием R. Вам следует попытаться получить исходный код R, создавший файл, вместо последнего измененного слоя.