Coverage for brodata / webservices.py: 65%

159 statements  

« prev     ^ index     » next       coverage.py v7.13.5, created at 2026-03-20 14:37 +0000

1import logging 

2import json 

3 

4import geopandas as gpd 

5import numpy as np 

6import requests 

7from requests.exceptions import HTTPError 

8from shapely.geometry import MultiPolygon, Point, Polygon 

9from .util import tqdm 

10 

11logger = logging.getLogger(__name__) 

12 

13 

14def get_gdf(kind, extent=None, config=None, index="DINO_NR", **kwargs): 

15 """Retrieve a GeoDataFrame for a configured BRO/GDN service. 

16 

17 Parameters 

18 ---------- 

19 kind : str 

20 Name of the service as defined in the configuration. 

21 extent : list, tuple, np.array, shapely.geometry.Polygon or shapely.geometry.MultiPolygon, optional 

22 Spatial filter. Use [xmin, xmax, ymin, ymax] for an envelope query or 

23 a polygon geometry for an intersects query. 

24 config : dict, optional 

25 Service configuration dictionary. If omitted, the default configuration 

26 from get_configuration() is used. 

27 index : str, optional 

28 Column name to use as index when present in the returned data. 

29 **kwargs 

30 Extra keyword arguments passed to requests.get via arcrest. 

31 

32 Returns 

33 ------- 

34 geopandas.GeoDataFrame 

35 GeoDataFrame with optional post-filtering and index assignment. 

36 """ 

37 if config is None: 

38 config = get_configuration() 

39 if kind not in config: 

40 raise (ValueError(f"Unknown kind: {kind}")) 

41 url = config[kind]["mapserver"] 

42 layer = 0 

43 if "layer" in config[kind]: 

44 layer = config[kind]["layer"] 

45 gdf = arcrest(url, layer, extent=extent, **kwargs) 

46 if "greater_than_0" in config[kind]: 

47 greater_than_0 = config[kind]["greater_than_0"] 

48 if greater_than_0 in gdf.columns: 

49 gdf = gdf[gdf[greater_than_0] > 0] 

50 if not gdf.empty and index in gdf.columns: 

51 gdf = gdf.set_index(index) 

52 return gdf 

53 

54 

55def arcrest( 

56 url, 

57 layer, 

58 extent=None, 

59 sr=28992, 

60 f="geojson", 

61 max_record_count=None, 

62 timeout=120, 

63 **kwargs, 

64): 

65 """Download data from an arcgis rest FeatureServer. 

66 

67 Parameters 

68 ---------- 

69 url : str 

70 arcrest url. 

71 layer : str 

72 layer 

73 extent : list, tuple, np.array, shapely.geometry.Polygon or shapely.geometry.MultiPolygon 

74 spatial filter as envelope [xmin, xmax, ymin, ymax] or polygon geometry. 

75 sr : int, optional 

76 coördinate reference system. The default is 28992 (RD). 

77 f : str, optional 

78 output format. Default is geojson 

79 max_record_count : int, optional 

80 maximum number of records for request. 

81 timeout : int, optional 

82 timeout time of request. Default is 120. 

83 

84 """ 

85 params = { 

86 "f": f, 

87 "outFields": "*", 

88 "outSR": sr, 

89 "where": "1=1", 

90 } 

91 if extent is not None: 

92 params["spatialRel"] = "esriSpatialRelIntersects" 

93 if isinstance(extent, (Polygon, MultiPolygon)): 

94 params["geometry"] = json.dumps(_as_arcgis_polygon(extent, sr=sr)) 

95 params["geometryType"] = "esriGeometryPolygon" 

96 else: 

97 xmin, xmax, ymin, ymax = extent 

98 params["geometry"] = f"{xmin},{ymin},{xmax},{ymax}" 

99 params["geometryType"] = "esriGeometryEnvelope" 

100 params["inSR"] = sr 

101 props = _get_data(url, {"f": "json"}, timeout=timeout, **kwargs) 

102 if max_record_count is None: 

103 max_record_count = props["maxRecordCount"] 

104 else: 

105 max_record_count = min(max_record_count, props["maxRecordCount"]) 

106 

107 params["returnIdsOnly"] = True 

108 url_query = f"{url}/{layer}/query" 

109 props = _get_data(url_query, params, timeout=timeout, **kwargs) 

110 params.pop("returnIdsOnly") 

111 if "objectIds" in props: 

112 object_ids = props["objectIds"] 

113 object_id_field_name = props["objectIdFieldName"] 

114 else: 

115 object_ids = props["properties"]["objectIds"] 

116 object_id_field_name = props["properties"]["objectIdFieldName"] 

117 if object_ids is not None and len(object_ids) > max_record_count: 

118 # download in batches 

119 object_ids.sort() 

120 n_d = int(np.ceil((len(object_ids) / max_record_count))) 

121 features = [] 

122 for i_d in tqdm(range(n_d)): 

123 i_min = i_d * max_record_count 

124 i_max = min(i_min + max_record_count - 1, len(object_ids) - 1) 

125 where = "{}>={} and {}<={}".format( 

126 object_id_field_name, 

127 object_ids[i_min], 

128 object_id_field_name, 

129 object_ids[i_max], 

130 ) 

131 params["where"] = where 

132 data = _get_data(url_query, params, timeout=timeout, **kwargs) 

133 features.extend(data["features"]) 

134 else: 

135 # download all data in one go 

136 data = _get_data(url_query, params, timeout=timeout, **kwargs) 

137 features = data["features"] 

138 if f == "json" or f == "pjson": 

139 # Interpret the geometry field 

140 data = [] 

141 for feature in features: 

142 if "rings" in feature["geometry"]: 

143 if len(feature["geometry"]) > 1: 

144 raise (NotImplementedError("Multiple rings not supported yet")) 

145 if len(feature["geometry"]["rings"]) == 1: 

146 geometry = Polygon(feature["geometry"]["rings"][0]) 

147 else: 

148 pols = [Polygon(xy) for xy in feature["geometry"]["rings"]] 

149 keep = [0] 

150 for i in range(1, len(pols)): 

151 if pols[i].within(pols[keep[-1]]): 

152 pols[keep[-1]] = pols[keep[-1]].difference(pols[i]) 

153 else: 

154 keep.append(i) 

155 if len(keep) == 1: 

156 geometry = pols[keep[0]] 

157 else: 

158 geometry = MultiPolygon([pols[i] for i in keep]) 

159 elif ( 

160 len(feature["geometry"]) == 2 

161 and "x" in feature["geometry"] 

162 and "y" in feature["geometry"] 

163 ): 

164 geometry = Point(feature["geometry"]["x"], feature["geometry"]["y"]) 

165 else: 

166 raise (Exception("Not supported yet")) 

167 feature["attributes"]["geometry"] = geometry 

168 data.append(feature["attributes"]) 

169 if len(data) == 0: 

170 # Assigning CRS to a GeoDataFrame without a geometry column is not supported 

171 gdf = gpd.GeoDataFrame() 

172 else: 

173 gdf = gpd.GeoDataFrame(data, crs=sr) 

174 else: 

175 # for geojson-data we can transform to GeoDataFrame right away 

176 if len(features) == 0: 

177 # Assigning CRS to a GeoDataFrame without a geometry column is not supported 

178 gdf = gpd.GeoDataFrame() 

179 else: 

180 gdf = gpd.GeoDataFrame.from_features(features, crs=sr) 

181 

182 return gdf 

183 

184 

185def _as_arcgis_polygon(polygon, sr=28992): 

186 """Convert polygon input to ArcGIS REST polygon geometry.""" 

187 if isinstance(polygon, dict): 

188 geometry = polygon.copy() 

189 if "rings" not in geometry: 

190 raise ValueError("polygon dict must contain 'rings'") 

191 elif isinstance(polygon, Polygon): 

192 rings = [list(polygon.exterior.coords)] 

193 for interior in polygon.interiors: 

194 rings.append(list(interior.coords)) 

195 geometry = {"rings": rings} 

196 elif isinstance(polygon, MultiPolygon): 

197 rings = [] 

198 for pol in polygon.geoms: 

199 rings.append(list(pol.exterior.coords)) 

200 for interior in pol.interiors: 

201 rings.append(list(interior.coords)) 

202 geometry = {"rings": rings} 

203 else: 

204 raise TypeError( 

205 "polygon must be a shapely Polygon, MultiPolygon, or dict with 'rings'" 

206 ) 

207 

208 if "spatialReference" not in geometry: 

209 geometry["spatialReference"] = {"wkid": sr} 

210 return geometry 

211 

212 

213def _get_data(url, params, timeout=5, **kwargs): 

214 """get data using a request 

215 

216 Parameters 

217 ---------- 

218 url : str 

219 url 

220 params : dict 

221 request parameters 

222 timeout : int, optional 

223 timeout time of request. Default is 120. 

224 

225 Returns 

226 ------- 

227 data 

228 

229 """ 

230 r = requests.get(url, params=params, timeout=timeout, **kwargs) 

231 if not r.ok: 

232 raise (HTTPError(f"Request not successful: {r.url}")) 

233 data = r.json() 

234 if "error" in data: 

235 code = data["error"]["code"] 

236 message = data["error"]["message"] 

237 raise (Exception(f"Error code {code}: {message}")) 

238 return data 

239 

240 

241services = { 

242 "lks_abo_rd": "Archeologisch booronderzoek", 

243 # "lks_aeo_rd": "", 

244 # "lks_aep_rd": "", 

245 "lks_bhr_g_rd": "Geologisch booronderzoek (BRO)", 

246 "lks_bhr_gt_rd": "Geotechnisch booronderzoek (BRO)", 

247 "lks_bhr_rd": "Bodemkundig booronderzoek (BRO)", 

248 "lks_cpt_rd": "Geotechnisch sondeeronderzoek (BRO)", 

249 # "lks_dgm_rd": "", 

250 # "lks_dgmdiep_rd": "", 

251 "lks_gbo_rd": "Geologisch booronderzoek", 

252 # "lks_geok_line_rd": "", 

253 # "lks_geok_rd": "", 

254 # "lks_geok_urania_rd": "", 

255 # "lks_gmd5_rd": "", 

256 # "lks_gmm_rd": "", 

257 "lks_gmw_rd": "Grondwatermonitoringput", 

258 "lks_gso_rd": "Geotechnisch sondeeronderzoek (GDN)", 

259 # "lks_gtp_rd": "", 

260 "lks_guf_rd": "Grondwatergebruiksystemen", 

261 "lks_gwo_rd": "Put met onderzoekgegevens", 

262 # "lks_mdl_rd": "", 

263 # "lks_nzs_rd": "", 

264 "lks_owo_rd": "Oppervlaktewateronderzoek", 

265 # "lks_rgs_rd": "", 

266 "lks_sfr_rd": "Bodemkundig wandonderzoek (BRO)", 

267 # "lks_sgm2_rd": "", 

268 "lks_vso_rd": "Verticaal elektrisch sondeeronderzoek", 

269 "lks_wbo_rd": "Geologisch waterbodemonderzoek", 

270 # "lks_wdm_brh_rd": "", 

271 # "lks_wdm_ghg_rd": "", 

272 # "lks_wdm_glg_rd": "", 

273 # "lks_wdm_gvg_rd": "", 

274 #'lks_wdm_gxg_rd": "", 

275 #'lks_wdm_ref_rd": "" 

276} 

277 

278 

279def get_configuration(mapserver_url=None): 

280 config = {} 

281 if mapserver_url is None: 

282 mapserver_url = "https://www.broloket.nl/standalone/rest/services" 

283 

284 # BRO 

285 bro_mapserver_url = f"{mapserver_url}/uitgifteloket_bro" 

286 config["Geologisch booronderzoek (BRO)"] = { 

287 "mapserver": f"{bro_mapserver_url}/lks_bhr_g_rd_v1/MapServer", 

288 "abbr": "bhrg", 

289 "rest_url": "https://publiek.broservices.nl/sr/bhrg/v3", 

290 } 

291 config["Geotechnisch booronderzoek (BRO)"] = { 

292 "mapserver": f"{bro_mapserver_url}/lks_bhr_gt_rd_v1/MapServer", 

293 "abbr": "bhrgt", 

294 "object": "BHR-GT", 

295 "rest_url": "https://publiek.broservices.nl/sr/bhrgt/v2", 

296 } 

297 config["Bodemkundig booronderzoek (BRO)"] = { 

298 "mapserver": f"{bro_mapserver_url}/lks_bhr_rd_v1/MapServer", 

299 "abbr": "bhrp", 

300 "object": "BHR_O", 

301 "rest_url": "https://publiek.broservices.nl/sr/bhrp/v2", 

302 } 

303 config["Geotechnisch sondeeronderzoek (BRO)"] = { 

304 "mapserver": f"{bro_mapserver_url}/lks_cpt_rd_v1/MapServer", 

305 "abbr": "cpt", 

306 "rest_url": "https://publiek.broservices.nl/sr/cpt/v1", 

307 } 

308 config["Bodemkundig wandonderzoek (BRO)"] = { 

309 "mapserver": f"{bro_mapserver_url}/lks_sfr_rd_v1/MapServer", 

310 "name_en": "pedological SoilFaceResearch", 

311 "abbr": "sfr", 

312 } 

313 config["Grondwatermonitoringput"] = { 

314 "mapserver": f"{bro_mapserver_url}/lks_gmw_rd_v1/MapServer", 

315 "rest_url": "https://publiek.broservices.nl/gm/gmw/v1", 

316 } 

317 config["Grondwaterstandonderzoek"] = { 

318 "reeks": "https://publiek.broservices.nl/gm/gld/v1/seriesAsCsv", 

319 "rest_url": "https://publiek.broservices.nl/gm/gld/v1", 

320 "abbr": "gld", 

321 } 

322 config["Grondwatersamenstellingsonderzoek"] = { 

323 "abbr": "gar", 

324 "rest_url": "https://publiek.broservices.nl/gm/gar/v1", 

325 "object": "GAR_O", 

326 } 

327 config["Grondwatergebruiksystemen"] = { 

328 "mapserver": f"{bro_mapserver_url}/lks_guf_rd_v1/MapServer", 

329 "name_nl": "Grondwatergebruiksystemen", 

330 "name_en": "Groundwater utilisation facility", 

331 "class": "GroundwaterUtilisationFacility", 

332 "abbr": "guf", 

333 "object": "GUF_PO", 

334 "rest_url": "https://publiek.broservices.nl/gu/guf/v1", 

335 } 

336 

337 # DINO 

338 gdn_mapserver_url = f"{mapserver_url}/uitgifteloket_gdn" 

339 config["Geologisch booronderzoek"] = { 

340 "mapserver": f"{gdn_mapserver_url}/lks_gbo_rd_v1/MapServer", 

341 "download": "https://www.dinoloket.nl/uitgifteloket/api/brh/sampledescription/csv", 

342 } 

343 config["Boormonsterprofiel"] = config["Geologisch booronderzoek"].copy() 

344 config["Boormonsterprofiel"]["greater_than_0"] = "MP_CNT" 

345 

346 config["Boormonsterfoto"] = config["Geologisch booronderzoek"].copy() 

347 config["Boormonsterfoto"]["greater_than_0"] = "MF_CNT" 

348 

349 config["Boorgatmeting"] = config["Geologisch booronderzoek"].copy() 

350 config["Boorgatmeting"]["greater_than_0"] = "BM_CNT" 

351 

352 config["Chemische analyse"] = config["Geologisch booronderzoek"].copy() 

353 config["Chemische analyse"]["greater_than_0"] = "CA_CNT" 

354 

355 config["Korrelgrootte analyse"] = config["Geologisch booronderzoek"].copy() 

356 config["Korrelgrootte analyse"]["greater_than_0"] = "KA_CNT" 

357 

358 config["Geologisch waterbodemonderzoek"] = { 

359 "mapserver": f"{gdn_mapserver_url}/lks_wbo_rd_v1/MapServer", 

360 "download": config["Geologisch booronderzoek"]["download"], 

361 } 

362 config["Archeologisch booronderzoek"] = { 

363 "mapserver": f"{gdn_mapserver_url}/lks_abo_rd_v1/MapServer", 

364 "download": config["Geologisch booronderzoek"]["download"], 

365 } 

366 config["Geotechnisch sondeeronderzoek (GDN)"] = { 

367 "mapserver": f"{gdn_mapserver_url}/lks_gso_rd_v1/MapServer", 

368 } 

369 config["Put met onderzoekgegevens"] = { 

370 "mapserver": f"{gdn_mapserver_url}/lks_gwo_rd_tiled_v1/MapServer", 

371 } 

372 config["Grondwatersamenstelling"] = { 

373 "mapserver": config["Put met onderzoekgegevens"]["mapserver"], 

374 "table": 2, # Grondwatersamenstelling 

375 "greater_than_0": "SA_CNT", 

376 "download": "https://www.dinoloket.nl/uitgifteloket/api/wo/gwo/qua/report", 

377 } 

378 config["Grondwaterstand"] = { 

379 "mapserver": config["Put met onderzoekgegevens"]["mapserver"], 

380 "table": 3, # Grondwaterstand 

381 "greater_than_0": "ST_CNT", 

382 "download": "https://www.dinoloket.nl/uitgifteloket/api/wo/gwo/full", 

383 "details": "https://www.dinoloket.nl/uitgifteloket/api/wo/gwo/details", 

384 } 

385 config["Oppervlaktewateronderzoek"] = { 

386 "mapserver": f"{gdn_mapserver_url}/lks_owo_rd_v1/MapServer" 

387 } 

388 config["Verticaal elektrisch sondeeronderzoek"] = { 

389 "mapserver": f"{gdn_mapserver_url}/lks_vso_rd_v1/MapServer", 

390 "download": "https://www.dinoloket.nl/uitgifteloket/api/ves/csv", 

391 } 

392 

393 return config