Добавить отсутствующий CRS в объект sf в R простых функциях - PullRequest
0 голосов
/ 20 апреля 2020

Есть ли способ добавить отсутствующий CRS? Мне дали несколько старых шейп-файлов Esri, и большая часть из них имеет CRS: NA, я хотел бы иметь возможность пакетной установки отсутствующего CRS. Я могу установить его вручную, когда импортирую в QGIS, но не могу установить его с помощью R. Я попытался st_transform(), но, похоже, он не работает, если в первую очередь не установлен CRS.

1 Ответ

0 голосов
/ 21 апреля 2020

Спасибо за комментарии. Я пробовал st_crs <- и st_set_crs(), но все еще не мог установить его правильно или построить шейп-файлы. Я заметил, что те, у которых отсутствовали .prj файлы, не имели CRS, что имеет смысл. поэтому я написал кусок для создания файла proj для каждого шейп-файла в каталоге, если он еще не существует. Немного хак, но это работает:)

РЕДАКТИРОВАТЬ: Добавлена ​​грубая система ведения журнала для отслеживания изменений

p = "/path/to/shape/files"

proj_str_bng = 'PROJCS["British_National_Grid",GEOGCS["GCS_OSGB_1936",DATUM["D_OSGB_1936",SPHEROID["Airy_1830",6377563.396,299.3249646]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Transverse_Mercator"],PARAMETER["False_Easting",400000.0],PARAMETER["False_Northing",-100000.0],PARAMETER["Central_Meridian",-2.0],PARAMETER["Scale_Factor",0.999601272],PARAMETER["Latitude_Of_Origin",49.0],UNIT["Meter",1.0]]'

proj_str_wgs84 = 'GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137.0,298.257223563]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]]'

shp_to_prj <- list.files(paste0(wd, p), pattern = "\\.shp$")

if (!file.exists(paste0(wd,p,"prj_log"))){
write("file, dt", file = paste0(wd,p,"prj_log"))
}

for (i in 1:length(shp_to_prj)) {
  shp_to_prj[[i]] <- paste0(wd, p, str_remove_all(shp_to_prj[[i]], ".shp"), ".prj")
  if(file.exists(shp_to_prj[[i]])) {
    next
  }
  write(proj_str_bng, file = shp_to_prj[[i]])
  write(paste0( shp_to_prj[[i]], ",",Sys.time()), file = paste0(wd,p,"prj_log"), append=TRUE)
}
write(paste0( "^===============================================================^", ",",Sys.time()), file = paste0(wd,p,"prj_log"), append=TRUE)

...