Code used to insert a polygon for Champaign County into a PostGIS geometry column.
Note that the sites.geometry column has three dimensions. Latitude and longitude are queried using ggplot2::map_data('counties') and elevation data is queried using rgbif::elevation(latidutde, longitude)
This query is used to define the boundary of site 1254, Champaign County, which will be used for testing and development of regional runs of PEcAn
library(ggplot2)
library(data.table)
library(rgbif)
##
## New to rgbif? Tutorial at http://ropensci.org/tutorials/rgbif_tutorial.html
## citation(package='rgbif') for the citation for this package
## Use suppressPackageStartupMessages() to suppress these startup messages in the future
d <- data.table(map_data("county"))
champaign <- d[region == "illinois" & subregion == "champaign"]
print(champaign)
## long lat group order region subregion
## 1: -88.46468 40.27893 570 22073 illinois champaign
## 2: -88.45895 40.38780 570 22074 illinois champaign
## 3: -88.01205 40.39352 570 22075 illinois champaign
## 4: -88.00632 40.37634 570 22076 illinois champaign
## 5: -87.93183 40.37061 570 22077 illinois champaign
## 6: -87.93183 40.22164 570 22078 illinois champaign
## 7: -87.94329 40.22164 570 22079 illinois champaign
## 8: -87.94329 39.87214 570 22080 illinois champaign
## 9: -88.47041 39.87214 570 22081 illinois champaign
## 10: -88.47614 40.21018 570 22082 illinois champaign
## 11: -88.46468 40.21591 570 22083 illinois champaign
## 12: -88.46468 40.27893 570 22084 illinois champaign
### for PostGIS, need polygon to return to start point
champaign <- rbind(champaign, champaign[1])
champaign_pts <- data.table(champaign[,elevation(latitude = lat, longitude = long) ])
print(champaign_pts)
## latitude longitude elevation
## 1: 40.27893 -88.46468 231.2208
## 2: 40.38780 -88.45895 225.4893
## 3: 40.39352 -88.01205 226.0480
## 4: 40.37634 -88.00632 227.0900
## 5: 40.37061 -87.93183 222.3453
## 6: 40.22164 -87.93183 224.1132
## 7: 40.22164 -87.94329 221.8404
## 8: 39.87214 -87.94329 218.7662
## 9: 39.87214 -88.47041 204.1542
## 10: 40.21018 -88.47614 223.6267
## 11: 40.21591 -88.46468 235.4709
## 12: 40.27893 -88.46468 231.2208
## 13: 40.27893 -88.46468 231.2208
boundary <- champaign_pts[,paste(longitude, latitude, elevation, collapse = ",")]
print(boundary)
## [1] "-88.4646835327148 40.2789344787598 231.220825195312,-88.4589538574219 40.3877983093262 225.489318847656,-88.0120468139648 40.3935241699219 226.048049926758,-88.0063171386719 40.376335144043 227.089981079102,-87.9318313598633 40.37060546875 222.345275878906,-87.9318313598633 40.2216377258301 224.113204956055,-87.9432907104492 40.2216377258301 221.840393066406,-87.9432907104492 39.8721351623535 218.766220092773,-88.4704132080078 39.8721351623535 204.154159545898,-88.4761428833008 40.2101783752441 223.626693725586,-88.4646835327148 40.2159080505371 235.470901489258,-88.4646835327148 40.2789344787598 231.220825195312,-88.4646835327148 40.2789344787598 231.220825195312"
#geometry <- paste0("ST_SetSRID(ST_MakePolygon(ST_Force_3D(ST_GeomFromText('LINESTRING(", boundary , ")'))), 4326)")
geometry <- paste0("ST_SetSRID(ST_MakePolygon(ST_GeomFromText('LINESTRING(", boundary , ")')), 4326)")
print(geometry)
## [1] "ST_SetSRID(ST_MakePolygon(ST_GeomFromText('LINESTRING(-88.4646835327148 40.2789344787598 231.220825195312,-88.4589538574219 40.3877983093262 225.489318847656,-88.0120468139648 40.3935241699219 226.048049926758,-88.0063171386719 40.376335144043 227.089981079102,-87.9318313598633 40.37060546875 222.345275878906,-87.9318313598633 40.2216377258301 224.113204956055,-87.9432907104492 40.2216377258301 221.840393066406,-87.9432907104492 39.8721351623535 218.766220092773,-88.4704132080078 39.8721351623535 204.154159545898,-88.4761428833008 40.2101783752441 223.626693725586,-88.4646835327148 40.2159080505371 235.470901489258,-88.4646835327148 40.2789344787598 231.220825195312,-88.4646835327148 40.2789344787598 231.220825195312)')), 4326)"
# insert statement
insert_polygon <- paste0("INSERT into sites (sitename, geometry) VALUES ('Champaign County', ", geometry, ");")
update_polygon <- paste0("update sites set geometry = ", geometry , " where sitename = 'Champaign County';")
print(insert_polygon)
## [1] "INSERT into sites (sitename, geometry) VALUES ('Champaign County', ST_SetSRID(ST_MakePolygon(ST_GeomFromText('LINESTRING(-88.4646835327148 40.2789344787598 231.220825195312,-88.4589538574219 40.3877983093262 225.489318847656,-88.0120468139648 40.3935241699219 226.048049926758,-88.0063171386719 40.376335144043 227.089981079102,-87.9318313598633 40.37060546875 222.345275878906,-87.9318313598633 40.2216377258301 224.113204956055,-87.9432907104492 40.2216377258301 221.840393066406,-87.9432907104492 39.8721351623535 218.766220092773,-88.4704132080078 39.8721351623535 204.154159545898,-88.4761428833008 40.2101783752441 223.626693725586,-88.4646835327148 40.2159080505371 235.470901489258,-88.4646835327148 40.2789344787598 231.220825195312,-88.4646835327148 40.2789344787598 231.220825195312)')), 4326));"
library(PEcAn.BETYdb)
settings <- read.settings('settings.xml')
db.query(insert_polygon, con = db.open(settings$database$bety))