Блог

У меня возникла проблема с сохранением растра, так как я получаю предупреждающее сообщение In .gd_SetProject(object, ...) : NOT UPDATED FOR PROJ >= 6

Я искал в Интернете и обновил все соответствующие пакеты из исходного кода и думаю, что это проблема с crs, однако я довольно новичок в R и не могу найти решение.

Вот мой код:

 # Create a template raster with cells equal in size to DEM
res(DEM)
aggregate <- raster::aggregate
r <- aggregate(DEM, 6)
res(r)

r <- SA %>%
st_transform(crs = 4326)%>%
rasterize(r, field = 1) %>%
# remove any empty cells
trim()

# Save raster:
r <- writeRaster(r, filename = "prediction-surface_SA.tif", overwrite = TRUE)
Warning message:
In .gd_SetProject(object, ...) : NOT UPDATED FOR PROJ >= 6

> DEM
class      : RasterLayer 
dimensions : 324, 435, 140940  (nrow, ncol, ncell)
resolution : 0.0008333333, 0.0008333333  (x, y)
extent     : 47.66083, 48.02333, -19.80917, -19.53917  (xmin, xmax, ymin, ymax)
crs        : NA 
source     : memory
names      : layer 
values     : 548, 1780  (min, max)

SA<-structure(list(FID = 0, geometry = structure(list(structure(list(
    structure(c(47.660604941588, 47.660604941588, 48.0232823584119, 
    48.0232823584119, 47.660604941588, -19.8088030171198, -19.5393336630522, 
    -19.5393336630522, -19.8088030171198, -19.8088030171198), .Dim = c(5L, 
    2L))), class = c("XY", "POLYGON", "sfg"))), n_empty = 0L, crs = structure(list(
    input = "WGS 84", wkt = "GEOGCRS["WGS 84",n    DATUM["World Geodetic System 1984",n        ELLIPSOID["WGS 84",6378137,298.257223563,n            LENGTHUNIT["metre",1]]],n    PRIMEM["Greenwich",0,n        ANGLEUNIT["degree",0.0174532925199433]],n    CS[ellipsoidal,2],n        AXIS["latitude",north,n            ORDER[1],n            ANGLEUNIT["degree",0.0174532925199433]],n        AXIS["longitude",east,n            ORDER[2],n            ANGLEUNIT["degree",0.0174532925199433]],n    ID["EPSG",4326]]"), class = "crs"), class = c("sfc_POLYGON", 
"sfc"), precision = 0, bbox = structure(c(xmin = 47.660604941588, 
ymin = -19.8088030171198, xmax = 48.0232823584119, ymax = -19.5393336630522
), class = "bbox"))), row.names = c(NA, -1L), class = c("sf", 
"data.frame"), sf_column = "geometry", agr = structure(c(FID = NA_integer_), class = "factor", .Label = c("constant", 
"aggregate", "identity")))

r<-new("RasterLayer", file = new(".RasterFile", name = "C:\Users\e77-smith\OneDrive - UWE Bristol (Staff)\Occupancy Models\Madagascar\prediction-surface_SA.tif", 
    datanotation = "FLT4S", byteorder = "little", nodatavalue = -Inf, 
    NAchanged = FALSE, nbands = 1L, bandorder = "BIL", offset = 0L, 
    toptobottom = TRUE, blockrows = 28L, blockcols = 72L, driver = "gdal", 
    open = FALSE), data = new(".SingleLayerData", values = logical(0), 
    offset = 0, gain = 1, inmemory = FALSE, fromdisk = TRUE, 
    isfactor = FALSE, attributes = list(), haveminmax = TRUE, 
    min = 1, max = 1, band = 1L, unit = "", names = "prediction.surface_SA"), 
    legend = new(".RasterLegend", type = character(0), values = logical(0), 
        color = logical(0), names = logical(0), colortable = logical(0)), 
    title = character(0), extent = new("Extent", xmin = 47.660833333, 
        xmax = 48.020833333, ymin = -19.809166667, ymax = -19.539166667), 
    rotated = FALSE, rotation = new(".Rotation", geotrans = numeric(0), 
        transfun = function () 
        NULL), ncols = 72L, nrows = 54L, crs = new("CRS", projargs = NA_character_), 
    history = list(), z = list())
 

Я думаю, что вы получаете это предупреждение, потому что система координат для DEM (и, следовательно r ,) не задана.

 library(raster)
r <- raster(res=30)
crs(r) <- NA
values(r) <- 1:ncell(r)
writeRaster(r, "tst.tif", overwrite=T)
#Warning message:
#In .gd_SetProject(object, ...) : NOT UPDATED FOR PROJ >= 6
 

Вы можете избежать предупреждения, установив crs

 crs(r) <- " proj=longlat"
writeRaster(r, "tst.tif", overwrite=T)