Return to Question

Corrected a code error that prevented the solution from running on my computer

Source Link

edited Dec 13, 2022 at 13:14

edit: setting mode = "wb" to download the file is crucial for correct reading using solution by Spacedman

I'm trying to get raster data on depth from EMODnet into R to analyse, using the WCS service.

library(RCurl)
library(XML)
library(raster)
xml <- paste(readLines("Request2.xml"),collapse = "\n")

myheader=c(Connection="close", 
            'Content-Type' = "application/xml",
             'Content-length'=nchar(xml))

response <- getURL(url = "http://ows.emodnet-bathymetry.eu/wcs?",
                   postfields = xml,
                   httpheader = myheader)

xmlresp <- xmlTreeParse(response, asText = TRUE,useInternalNodes=TRUE)
link <- xmlGetAttr(getNodeSet(xmlresp,'//wcs:Coverages/wcs:Coverage/ows:Reference' )[[1]],
           "href")
download.file(link, destfile = "test.tiff",
              mode = "wb")

ras <- terra::rast("test.tiff")

I'm trying to get raster data on depth from EMODnet into R to analyse, using the WCS service.

library(RCurl)
library(XML)
library(raster)
xml <- paste(readLines("Request2.xml"),collapse = "\n")

myheader=c(Connection="close", 
            'Content-Type' = "application/xml",
             'Content-length'=nchar(xml))

response <- getURL(url = "http://ows.emodnet-bathymetry.eu/wcs?",
                   postfields = xml,
                   httpheader = myheader)

xmlresp <- xmlTreeParse(response, asText = TRUE,useInternalNodes=TRUE)
link <- xmlGetAttr(getNodeSet(xmlresp,'//wcs:Coverages/wcs:Coverage/ows:Reference' )[[1]],
           "href")
download.file(link, destfile = "test.tiff")

ras <- terra::rast("test.tiff")

edit: setting mode = "wb" to download the file is crucial for correct reading using solution by Spacedman

I'm trying to get raster data on depth from EMODnet into R to analyse, using the WCS service.

library(RCurl)
library(XML)
library(raster)
xml <- paste(readLines("Request2.xml"),collapse = "\n")

myheader=c(Connection="close", 
            'Content-Type' = "application/xml",
             'Content-length'=nchar(xml))

response <- getURL(url = "http://ows.emodnet-bathymetry.eu/wcs?",
                   postfields = xml,
                   httpheader = myheader)

xmlresp <- xmlTreeParse(response, asText = TRUE,useInternalNodes=TRUE)
link <- xmlGetAttr(getNodeSet(xmlresp,'//wcs:Coverages/wcs:Coverage/ows:Reference' )[[1]],
           "href")
download.file(link, destfile = "test.tiff",
              mode = "wb")

ras <- terra::rast("test.tiff")

Source Link

asked Dec 12, 2022 at 15:49

Joris Meys

Process WCS data from EMODnet using R

I'm trying to get raster data on depth from EMODnet into R to analyse, using the WCS service.

Following the examples here, I've constructed the following file called "request2.xml":

<?xml version="1.0" encoding="UTF-8"?>
<wcs:GetCoverage xmlns:wcs="http://www.opengis.net/wcs/1.1.1" xmlns:gml="http://www.opengis.net/gml"
  xmlns:ogc="http://www.opengis.net/ogc" xmlns:ows="http://www.opengis.net/ows/1.1"
  xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance"
  version="1.1.1" service="WCS"
  xsi:schemaLocation="http://www.opengis.net/wcs/1.1.1 http://schemas.opengis.net/wcs/1.1.1/wcsAll.xsd">
   <ows:Identifier>emodnet:mean</ows:Identifier>
   <wcs:DomainSubset>
      <ows:BoundingBox crs="urn:ogc:def:crs:EPSG::4326">
         <ows:LowerCorner>43.48729 -2.5319802324929</ows:LowerCorner>
         <ows:UpperCorner>44.016457 -1.634064</ows:UpperCorner>
      </ows:BoundingBox>
   </wcs:DomainSubset>
   <wcs:Output store="true" format="image/tiff">
      <wcs:GridCRS>
         <wcs:GridBaseCRS>urn:ogc:def:crs:EPSG::4326</wcs:GridBaseCRS>
         <wcs:GridType>urn:ogc:def:method:WCS:1.1:2dSimpleGrid</wcs:GridType>
         <wcs:GridOffsets>0.0020881772848672093 -0.0020833346456692954</wcs:GridOffsets>
         <wcs:GridCS>urn:ogc:def:cs:OGC:0.0:Grid2dSquareCS</wcs:GridCS>
      </wcs:GridCRS>
   </wcs:Output>
</wcs:GetCoverage>

I try to process this using the following code in R:

library(RCurl)
library(XML)
library(raster)
xml <- paste(readLines("Request2.xml"),collapse = "\n")

myheader=c(Connection="close", 
            'Content-Type' = "application/xml",
             'Content-length'=nchar(xml))

response <- getURL(url = "http://ows.emodnet-bathymetry.eu/wcs?",
                   postfields = xml,
                   httpheader = myheader)

xmlresp <- xmlTreeParse(response, asText = TRUE,useInternalNodes=TRUE)
link <- xmlGetAttr(getNodeSet(xmlresp,'//wcs:Coverages/wcs:Coverage/ows:Reference' )[[1]],
           "href")
download.file(link, destfile = "test.tiff")

ras <- terra::rast("test.tiff")

Everything works, until I try to extract the raster information. I keep getting a message that the data is rotated, and I should use rectify(). But doing so loses all information. Furthermore, upon inspection it seems that the information isn't read correctly:

> ras <- terra::rast("test.tiff")
Warning: [rast] the data in this file are rotated. Use 'rectify' to fix that> 
> ras
class       : SpatRaster 
dimensions  : 431, 253, 1  (nrow, ncol, nlyr)
resolution  : 0, 0  (x, y)
extent      : -1.634064, -1.634064, 43.48729, 43.48729  (xmin, xmax, ymin, ymax)
coord. ref. : lon/lat WGS 84 (EPSG:4326) 
source      : test.tiff 
name        : test

I've looked at drivers for gdal and feel a bit like I'm drowning. Anywone who can point out what I'm doing wrong?

raster r wcs