Haribo
Reputation: 2226
transform a polygon class sf to class owin in r
My polygon has the following structure:
t1 = structure(list(WindfarmId = "NO39", Name = "Sørlige Nordsjø I",
Country = "NO", Status = 99L, WindfarmStatus = "Development Zone",
StatusComments = NA_character_, CapacityMWMin = 500, CapacityMWMax = 1500,
NoTurbinesMin = NA_integer_, NoTurbinesMax = NA_integer_,
Comments = "In 2013, NVE divided considered offshore wind zones into three categories. This area is considered Category A.\r\n\r\nCategory A: Wind power development within the zone \r\nis technically and economically feasible, and will have \r\nrelatively few negative impacts. Grid connection is possible \r\nbefore 2025. \r\n\r\nCategory B: Wind power development within the zone will \r\nhave challenges related to either technical aspects or conflict of interests/negative impacts. The challenges might \r\nbe resolved in the future through technology development, \r\ngrid measures and/or mitigation measures. NVE considers \r\nthat zones in this category can be opened when technology \r\nmatures, or when existing use of the areas changes. \r\n\r\nCategory C: Wind power development within the zone represents greater challenges than in the other two categories. \r\nConflicts of interest in the areas are not easily resolved. \r\nForeseen negative impacts are still considered acceptable. Zones in this category should not be opened at the \r\nexpense of zones in the two other categories.\r\n\r\nA 2018 review of the areas put Utsira Nord, and Sørlige Nordsjø I/II forward as recommended areas for opening up for offshore wind. The 2019 consultation hearing however did not include Sørlige Nordsjø I. Although Sørlige Nordsjø I has less conflict with shipping and lower environmental risk, Sørlige Nordsjø II is closer to the European power system, is a larger area and has fewer allocated petroleum production licences.",
TurbineMWMin = NA_real_, TurbineMWMax = NA_real_, OtherNames = "Offshore wind - study area; Category A",
CountryName = "Norway", Lat = 57.4203547589749, Lon = 3.5335822891109,
IsEstimatedLocation = 0L, IsOnHold = 0L, OffshoreConsStartDate = structure(NA_real_, class = c("POSIXct",
"POSIXt"), tzone = ""), ProjFullCommDate = structure(NA_real_, class = c("POSIXct",
"POSIXt"), tzone = ""), Georegion = "Europe", Developers = "Norges vassdrags- og energidirektorat",
Owners = NA_character_, Operators = NA_character_, WaterDepthMinM = 58,
WaterDepthMaxM = 72, GISFILE = "EU_SCRE_WindFarms_PN_4C_230110_4258",
Shape_Length = 2.25250768130164, Shape_Area = 0.203134377892303,
Shape = structure(list(structure(list(list(structure(c(4.02249908400006,
3.14500045800003, 3.04527854900005, 3.92694473300003, 4.02249908400006,
57.3724994670001, 57.250833512, 57.470277787, 57.5883331310001,
57.3724994670001), .Dim = c(5L, 2L)))), class = c("XY", "MULTIPOLYGON",
"sfg"))), class = c("sfc_MULTIPOLYGON", "sfc"), precision = 0, bbox = structure(c(xmin = 3.04527854900005,
ymin = 57.250833512, xmax = 4.02249908400006, ymax = 57.5883331310001
), class = "bbox"), crs = structure(list(input = "ETRS89",
wkt = "GEOGCRS[\"ETRS89\",\n DATUM[\"European Terrestrial Reference System 1989\",\n ELLIPSOID[\"GRS 1980\",6378137,298.257222101,\n LENGTHUNIT[\"metre\",1]]],\n PRIMEM[\"Greenwich\",0,\n ANGLEUNIT[\"degree\",0.0174532925199433]],\n CS[ellipsoidal,2],\n AXIS[\"geodetic latitude (Lat)\",north,\n ORDER[1],\n ANGLEUNIT[\"degree\",0.0174532925199433]],\n AXIS[\"geodetic longitude (Lon)\",east,\n ORDER[2],\n ANGLEUNIT[\"degree\",0.0174532925199433]],\n USAGE[\n SCOPE[\"unknown\"],\n AREA[\"Europe - ETRS89\"],\n BBOX[32.88,-16.1,84.17,40.18]],\n ID[\"EPSG\",4258]]"), class = "crs"), n_empty = 0L),
label = list(structure("Name: Sørlige Nordsjø I<br>Status: Development Zone<br>Capacity: 500[MW]<br>Owners: NA<br>", html = TRUE, class = c("html",
"character")))), sf_column = "Shape", agr = structure(c(WindfarmId = NA_integer_,
Name = NA_integer_, Country = NA_integer_, Status = NA_integer_,
WindfarmStatus = NA_integer_, StatusComments = NA_integer_, CapacityMWMin = NA_integer_,
CapacityMWMax = NA_integer_, NoTurbinesMin = NA_integer_, NoTurbinesMax = NA_integer_,
Comments = NA_integer_, TurbineMWMin = NA_integer_, TurbineMWMax = NA_integer_,
OtherNames = NA_integer_, CountryName = NA_integer_, Lat = NA_integer_,
Lon = NA_integer_, IsEstimatedLocation = NA_integer_, IsOnHold = NA_integer_,
OffshoreConsStartDate = NA_integer_, ProjFullCommDate = NA_integer_,
Georegion = NA_integer_, Developers = NA_integer_, Owners = NA_integer_,
Operators = NA_integer_, WaterDepthMinM = NA_integer_, WaterDepthMaxM = NA_integer_,
GISFILE = NA_integer_, Shape_Length = NA_integer_, Shape_Area = NA_integer_,
label = NA_integer_), .Label = c("constant", "aggregate", "identity"
), class = "factor"), row.names = 1287L, class = c("sf", "data.frame"
))
I'm trying to use the spatstat.geom package to convert t1 to class owin object :
t1_owin = as.owin(t1)
Error in as.owin.sfc_MULTIPOLYGON(st_cast(W, "MULTIPOLYGON"), ...) :
Only projected coordinates may be converted to spatstat class objects
And the coordinates in t1 are lat/lng :
st_coordinates(t1)
X Y L1 L2 L3
[1,] 4.022499 57.37250 1 1 1
[2,] 3.145000 57.25083 1 1 1
[3,] 3.045279 57.47028 1 1 1
[4,] 3.926945 57.58833 1 1 1
[5,] 4.022499 57.37250 1 1 1
How do I project the t1 to UTM in order to be able to get my owin object ?!
Upvotes: 0
Views: 192
Answers (1)
dieghernan
Reputation: 3402
You can use crsuggest to get the right CRS (also the corresponding UTC if you need) and project your sf object so you can make your owin. Full workflow:
t1 = structure(list(WindfarmId = "NO39", Name = "Sørlige Nordsjø I",
Country = "NO", Status = 99L, WindfarmStatus = "Development Zone",
StatusComments = NA_character_, CapacityMWMin = 500, CapacityMWMax = 1500,
NoTurbinesMin = NA_integer_, NoTurbinesMax = NA_integer_,
Comments = "In 2013, NVE divided considered offshore wind zones into three categories. This area is considered Category A.\r\n\r\nCategory A: Wind power development within the zone \r\nis technically and economically feasible, and will have \r\nrelatively few negative impacts. Grid connection is possible \r\nbefore 2025. \r\n\r\nCategory B: Wind power development within the zone will \r\nhave challenges related to either technical aspects or conflict of interests/negative impacts. The challenges might \r\nbe resolved in the future through technology development, \r\ngrid measures and/or mitigation measures. NVE considers \r\nthat zones in this category can be opened when technology \r\nmatures, or when existing use of the areas changes. \r\n\r\nCategory C: Wind power development within the zone represents greater challenges than in the other two categories. \r\nConflicts of interest in the areas are not easily resolved. \r\nForeseen negative impacts are still considered acceptable. Zones in this category should not be opened at the \r\nexpense of zones in the two other categories.\r\n\r\nA 2018 review of the areas put Utsira Nord, and Sørlige Nordsjø I/II forward as recommended areas for opening up for offshore wind. The 2019 consultation hearing however did not include Sørlige Nordsjø I. Although Sørlige Nordsjø I has less conflict with shipping and lower environmental risk, Sørlige Nordsjø II is closer to the European power system, is a larger area and has fewer allocated petroleum production licences.",
TurbineMWMin = NA_real_, TurbineMWMax = NA_real_, OtherNames = "Offshore wind - study area; Category A",
CountryName = "Norway", Lat = 57.4203547589749, Lon = 3.5335822891109,
IsEstimatedLocation = 0L, IsOnHold = 0L, OffshoreConsStartDate = structure(NA_real_, class = c("POSIXct",
"POSIXt"), tzone = ""), ProjFullCommDate = structure(NA_real_, class = c("POSIXct",
"POSIXt"), tzone = ""), Georegion = "Europe", Developers = "Norges vassdrags- og energidirektorat",
Owners = NA_character_, Operators = NA_character_, WaterDepthMinM = 58,
WaterDepthMaxM = 72, GISFILE = "EU_SCRE_WindFarms_PN_4C_230110_4258",
Shape_Length = 2.25250768130164, Shape_Area = 0.203134377892303,
Shape = structure(list(structure(list(list(structure(c(4.02249908400006,
3.14500045800003, 3.04527854900005, 3.92694473300003, 4.02249908400006,
57.3724994670001, 57.250833512, 57.470277787, 57.5883331310001,
57.3724994670001), .Dim = c(5L, 2L)))), class = c("XY", "MULTIPOLYGON",
"sfg"))), class = c("sfc_MULTIPOLYGON", "sfc"), precision = 0, bbox = structure(c(xmin = 3.04527854900005,
ymin = 57.250833512, xmax = 4.02249908400006, ymax = 57.5883331310001
), class = "bbox"), crs = structure(list(input = "ETRS89",
wkt = "GEOGCRS[\"ETRS89\",\n DATUM[\"European Terrestrial Reference System 1989\",\n ELLIPSOID[\"GRS 1980\",6378137,298.257222101,\n LENGTHUNIT[\"metre\",1]]],\n PRIMEM[\"Greenwich\",0,\n ANGLEUNIT[\"degree\",0.0174532925199433]],\n CS[ellipsoidal,2],\n AXIS[\"geodetic latitude (Lat)\",north,\n ORDER[1],\n ANGLEUNIT[\"degree\",0.0174532925199433]],\n AXIS[\"geodetic longitude (Lon)\",east,\n ORDER[2],\n ANGLEUNIT[\"degree\",0.0174532925199433]],\n USAGE[\n SCOPE[\"unknown\"],\n AREA[\"Europe - ETRS89\"],\n BBOX[32.88,-16.1,84.17,40.18]],\n ID[\"EPSG\",4258]]"), class = "crs"), n_empty = 0L),
label = list(structure("Name: Sørlige Nordsjø I<br>Status: Development Zone<br>Capacity: 500[MW]<br>Owners: NA<br>", html = TRUE, class = c("html",
"character")))), sf_column = "Shape", agr = structure(c(WindfarmId = NA_integer_,
Name = NA_integer_, Country = NA_integer_, Status = NA_integer_,
WindfarmStatus = NA_integer_, StatusComments = NA_integer_, CapacityMWMin = NA_integer_,
CapacityMWMax = NA_integer_, NoTurbinesMin = NA_integer_, NoTurbinesMax = NA_integer_,
Comments = NA_integer_, TurbineMWMin = NA_integer_, TurbineMWMax = NA_integer_,
OtherNames = NA_integer_, CountryName = NA_integer_, Lat = NA_integer_,
Lon = NA_integer_, IsEstimatedLocation = NA_integer_, IsOnHold = NA_integer_,
OffshoreConsStartDate = NA_integer_, ProjFullCommDate = NA_integer_,
Georegion = NA_integer_, Developers = NA_integer_, Owners = NA_integer_,
Operators = NA_integer_, WaterDepthMinM = NA_integer_, WaterDepthMaxM = NA_integer_,
GISFILE = NA_integer_, Shape_Length = NA_integer_, Shape_Area = NA_integer_,
label = NA_integer_), .Label = c("constant", "aggregate", "identity"
), class = "factor"), row.names = 1287L, class = c("sf", "data.frame"
))
library(sf)
#> Warning: package 'sf' was built under R version 4.1.3
#> Linking to GEOS 3.10.2, GDAL 3.4.1, PROJ 7.2.1; sf_use_s2() is TRUE
library(spatstat.geom)
#> Warning: package 'spatstat.geom' was built under R version 4.1.3
#> Loading required package: spatstat.data
#> Warning: package 'spatstat.data' was built under R version 4.1.3
#> spatstat.geom 3.0-3
as.owin(t1)
#> Error in as.owin.sfc_MULTIPOLYGON(st_cast(W, "MULTIPOLYGON"), ...): Only projected coordinates may be converted to spatstat class objects
# This is because
st_is_longlat(t1)
#> [1] TRUE
# So use a projection
candidates <- crsuggest::suggest_crs(t1)
candidates
#> # A tibble: 10 x 6
#> crs_code crs_name crs_type crs_gcs crs_units crs_pr~1
#> <chr> <chr> <chr> <dbl> <chr> <chr>
#> 1 3043 ETRS89 / UTM zone 31N (N-E) projected 4258 m +proj=u~
#> 2 25831 ETRS89 / UTM zone 31N projected 4258 m +proj=u~
#> 3 23031 ED50 / UTM zone 31N projected 4230 m +proj=u~
#> 4 3035 ETRS89-extended / LAEA Europe projected 4258 m +proj=l~
#> 5 3034 ETRS89-extended / LCC Europe projected 4258 m +proj=l~
#> 6 32631 WGS 84 / UTM zone 31N projected 4326 m +proj=u~
#> 7 32431 WGS 72BE / UTM zone 31N projected 4324 m +proj=u~
#> 8 32231 WGS 72 / UTM zone 31N projected 4322 m +proj=u~
#> 9 3832 WGS 84 / PDC Mercator projected 4326 m +proj=m~
#> 10 3576 WGS 84 / North Pole LAEA Russia projected 4326 m +proj=l~
#> # ... with abbreviated variable name 1: crs_proj4
# For the sake of example, I just choose the top 1, but you may want another
t1_proj <- st_transform(t1, as.numeric(candidates$crs_code[1]))
t1_proj
#> Simple feature collection with 1 feature and 31 fields
#> Geometry type: MULTIPOLYGON
#> Dimension: XY
#> Bounding box: xmin: 502715.8 ymin: 6345317 xmax: 561491.4 ymax: 6383258
#> Projected CRS: ETRS89 / UTM zone 31N (N-E)
#> WindfarmId Name Country Status WindfarmStatus
#> 1287 NO39 Sørlige Nordsjø I NO 99 Development Zone
#> StatusComments CapacityMWMin CapacityMWMax NoTurbinesMin NoTurbinesMax
#> 1287 <NA> 500 1500 NA NA
#> Comments
#> 1287 In 2013, NVE divided considered offshore wind zones into three categories. This area is considered Category A.\r\n\r\nCategory A: Wind power development within the zone \r\nis technically and economically feasible, and will have \r\nrelatively few negative impacts. Grid connection is possible \r\nbefore 2025. \r\n\r\nCategory B: Wind power development within the zone will \r\nhave challenges related to either technical aspects or conflict of interests/negative impacts. The challenges might \r\nbe resolved in the future through technology development, \r\ngrid measures and/or mitigation measures. NVE considers \r\nthat zones in this category can be opened when technology \r\nmatures, or when existing use of the areas changes. \r\n\r\nCategory C: Wind power development within the zone represents greater challenges than in the other two categories. \r\nConflicts of interest in the areas are not easily resolved. \r\nForeseen negative impacts are still considered acceptable. Zones in this category should not be opened at the \r\nexpense of zones in the two other categories.\r\n\r\nA 2018 review of the areas put Utsira Nord, and Sørlige Nordsjø I/II forward as recommended areas for opening up for offshore wind. The 2019 consultation hearing however did not include Sørlige Nordsjø I. Although Sørlige Nordsjø I has less conflict with shipping and lower environmental risk, Sørlige Nordsjø II is closer to the European power system, is a larger area and has fewer allocated petroleum production licences.
#> TurbineMWMin TurbineMWMax OtherNames
#> 1287 NA NA Offshore wind - study area; Category A
#> CountryName Lat Lon IsEstimatedLocation IsOnHold
#> 1287 Norway 57.42035 3.533582 0 0
#> OffshoreConsStartDate ProjFullCommDate Georegion
#> 1287 <NA> <NA> Europe
#> Developers Owners Operators WaterDepthMinM
#> 1287 Norges vassdrags- og energidirektorat <NA> <NA> 58
#> WaterDepthMaxM GISFILE Shape_Length Shape_Area
#> 1287 72 EU_SCRE_WindFarms_PN_4C_230110_4258 2.252508 0.2031344
#> Shape
#> 1287 MULTIPOLYGON (((561491.4 63...
#> label
#> 1287 Name: Sørlige Nordsjø I<br>Status: Development Zone<br>Capacity: 500[MW]<br>Owners: NA<br>
# And now
as.owin(t1_proj)
#> window: polygonal boundary
#> enclosing rectangle: [502715.8, 561491.4] x [6345317, 6383258] units
Created on 2023-02-02 with reprex v2.0.2
Upvotes: 3
Related Questions
- PostgresSQL: Create a function to update a table value
- How to implement dynamic sql function in postgresql?
- Postgresql function invocation
- Update a column in a function using plpgsql
- Update a table column using postgres function
- POSTGRESQL Function to update a field in table
- Calling a function for each updated row in postgresql
- Creating Postgres/PostGIS function to update table values based on spatial query
- Update Statement Using An Existing Field and A Function
- Update function in PostgreSQL