I'm using the following Python code to download a tiff file from a WMS server:
import requests
import numpy as np
import rasterio
from rasterio.transform import from_origin
from rasterio.transform import rowcol
from rasterio.transform import rowcol
from rasterio.merge import merge
from rasterio.plot import show
from pyproj import Transformer
import certifi
import ssl
import matplotlib.pyplot as plt
Layer = "warme_nachten_huidig"
minx=13588.806396509579
miny=306893.81943950325
maxx=277788.8063965096
maxy=612693.8194395032
map_width = "2642"
map_height = "3058"
bbox = (minx, miny, maxx, maxy)
def map_downloaden(bbox, width, height):
# WMS Service URL and parameters
ssl.SSLContext.verify_mode = ssl.VerifyMode.CERT_OPTIONAL
# Define the WMS service URL
wms_url = "?"
# Define the parameters for the GetMap request
params = {
"service": "WMS",
"version": "1.3.0",
"request": "GetMap",
"layers": f"{Layer}",
"styles": "",
"crs": "EPSG:28992",
"bbox": ','.join(map(str, bbox)), # Bounding box in EPSG:28992
"width": width, # Width of the output image
"height": height, # Height of the output image
"format": "image/geotiff" # Requesting GeoTIFF format
}
# Make the request to the WMS server
response = requests.get(wms_url, params=params, stream=True, verify=certifi.where())
# Check if the request was successful
if response.status_code == 200:
with open(file_path_python, 'wb') as f:
f.write(response.content)
print(f"Map successfully downloaded to {file_path_python}")
else:
print(f"Failed to download map. Status code: {response.status_code}")
map_downloaden(bbox, map_width, map_height)
When I compare the value of a specific coordinate of this file with the original file by using the following code:
def get_coordinate_value(coord, file_path, versie):
with rasterio.open(file_path) as dataset:
# Transform geographic coordinates to the raster's CRS
transformer = Transformer.from_crs("EPSG:4326", dataset.crs, always_xy=True)
projected_coord = transformer.transform(coord[0], coord[1])
# Transform geographic coordinates to raster's pixel coordinates
row, col = rowcol(dataset.transform, projected_coord[0], projected_coord[1])
# Ensure the calculated row and col are within the raster's dimensions
if 0 <= row < dataset.height and 0 <= col < dataset.width:
# Read the data value at the given pixel location
data_value = dataset.read(1)[row, col]
print(f"---- {versie}-------- Height: {dataset.height} Width: {dataset.width} ---- row: {row} col: {col} ---- COORDINATES: {coord} ---- VALUE: {data_value} ")
else:
print("Coordinates are out of the raster's bounds.")
coord = ( 5.975079577360298, 51.6954090242412)
get_coordinate_value(coord, file_path_python, "Python---")
get_coordinate_value(coord, file_path_original, "Original")
I get different values:
---- Python-----------
Height: 3058 Width: 2642 ---- row: 2006 col: 1820 ---- COORDINATES: (5.975079577360298, 51.6954090242412) ---- VALUE: 3
---- Original---------
Height: 3058 Width: 2642 ---- row: 2006 col: 1820 ---- COORDINATES: (5.975079577360298, 51.6954090242412) ---- VALUE: 0.6150805354118347
The original file is one that I downloaded from the website, the Python file is the one I am trying to retrieve from the WMS server.
Does anybody know why this happens? What am I doing wrong?