-
Notifications
You must be signed in to change notification settings - Fork 1
/
grb.py
545 lines (491 loc) · 20.1 KB
/
grb.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
import json
import logging
from copy import deepcopy
from datetime import date
from datetime import datetime
import numpy as np
from shapely import intersects
from shapely.geometry import shape
from brdr.aligner import Aligner
from brdr.constants import (
DEFAULT_CRS,
LAST_VERSION_DATE,
DATE_FORMAT,
VERSION_DATE,
FORMULA_FIELD_NAME,
BASE_FORMULA_FIELD_NAME,
)
from brdr.constants import DOWNLOAD_LIMIT
from brdr.constants import GRB_BUILDING_ID
from brdr.constants import GRB_FEATURE_URL
from brdr.constants import GRB_FISCAL_PARCELS_URL
from brdr.constants import GRB_KNW_ID
from brdr.constants import GRB_MAX_REFERENCE_BUFFER
from brdr.constants import GRB_PARCEL_ID
from brdr.constants import GRB_VERSION_DATE
from brdr.enums import GRBType, AlignerResultType
from brdr.geometry_utils import buffer_pos, safe_intersection, safe_unary_union
from brdr.geometry_utils import create_donut
from brdr.geometry_utils import features_by_geometric_operation
from brdr.geometry_utils import get_bbox
from brdr.loader import GeoJsonLoader, DictLoader
from brdr.logger import Logger
from brdr.utils import geojson_to_dicts
from brdr.utils import get_collection
from brdr.utils import get_collection_by_partition
log = logging.getLogger(__name__)
datetime_format_TZ = "%Y-%m-%dT%H:%M:%SZ"
def is_grb_changed(
geometry,
grb_type=GRBType.ADP,
date_start=date.today(),
date_end=date.today(),
border_distance=0,
crs=DEFAULT_CRS,
):
"""
checks if a geometry is possibly affected by changes in the reference layer
during a specified timespan
Args:
geometry: Geometry to check on GRB-changes
grb_type: Type of GRB (parcels, buildings,artwork,...) to check
date_start: Start of timespan to check if GRB-changes has occurred
date_end: End of timespan to check if GRB-changes has occurred
border_distance: Distance that can be used to only check the 'border' of
the geometry, so 'big' geometries with internal parcel-updates are not
affected
(Default:0, indicating that the full geometry is checked fot GRB-changes)
crs: Coordinate reference system to use
Returns: Boolean, indicating if GRB is changed underneath the geometry
"""
last_version_date = get_last_version_date(
geometry, grb_type=grb_type, crs=crs, border_distance=border_distance
)
if last_version_date is None:
return None
if date_start <= last_version_date <= date_end:
return True
return False
def get_affected_by_grb_change(
dict_thematic,
grb_type=GRBType.ADP,
date_start=date.today(),
date_end=date.today(),
one_by_one=False,
border_distance=0,
geometry_thematic_union=None,
crs=DEFAULT_CRS,
):
"""
Get a list of affected and unaffected IDs by GRB-changes in a
specific timespan
Args:
dict_thematic: dictionary if thematicID & Geometry
grb_type: Type of GRB: parcels, buildings,...
date_start: start-date to check changes in GRB
date_end: end-date to check changes in GRB
one_by_one: parameter to choose the methodology to check changes:
* True: Every thematic geometry is checked individually (loop)
* False: All GRB-parcels intersecting the thematic dictionary is checked
at once
border_distance: Distance that can be used to only check the 'border' of the
geometry, so 'big' geometries with internal parcel-updates are not affected
(Default:0, indicating that the full geometry is checked fot GRB-changes)
Returns:
dictionary of affected geometries
"""
affected = []
unaffected = []
if border_distance > 0:
for key in dict_thematic.keys():
dict_thematic[key] = create_donut(dict_thematic[key], border_distance)
if one_by_one:
for key in dict_thematic:
geom = dict_thematic[key]
if is_grb_changed(geom, grb_type, date_start, date_end):
affected.append(key)
else:
unaffected.append(key)
return affected, unaffected
else:
# Temporal filter on VERDATUM
if geometry_thematic_union is None:
geometry_thematic_union = safe_unary_union(list(dict_thematic.values()))
coll_changed_grb, name_reference_id = get_collection_grb_actual(
geometry_thematic_union,
grb_type=grb_type,
partition=1000,
date_start=date_start,
date_end=date_end,
crs=crs,
)
dict_changed_grb, dict_changed_grb_properties = geojson_to_dicts(
coll_changed_grb, name_reference_id
)
if len(dict_changed_grb) == 0:
logging.info("No detected changes")
return affected, list(dict_thematic.keys()) # empty affected dict
logging.info("Changed parcels in timespan: " + str(len(dict_changed_grb)))
thematic_intersections = features_by_geometric_operation(
list(dict_thematic.values()),
list(dict_thematic.keys()),
list(dict_changed_grb.values()),
predicate="intersects",
)
logging.info("Number of filtered features: " + str(len(thematic_intersections)))
for key, geom in dict_thematic.items():
(
affected.append(key)
if key in thematic_intersections
else unaffected.append(key)
)
return affected, unaffected
def get_last_version_date(
geometry,
grb_type=GRBType.ADP,
crs=DEFAULT_CRS,
limit=DOWNLOAD_LIMIT,
border_distance=0,
):
"""
Retrieves the date of the last version for a specific geographic area within
GRB (parcels, buildings,...).
This function queries the GRB-API to find the most recent version-date (=last
update of object) for reference data of the specified `grb_type` (e.g., ADP,
GBG, KNW) within the boundary of the provided `geometry`.
Args:
geometry (BaseGeometry): A Shapely geometry representing the area of interest.
grb_type (GRBType, optional): The type of GRB data to consider.
Defaults to GRBType.ADP (administrative plots).
Returns:
str: The date of the last version for the specified GRB data type within the area,
formatted as a string according to the GRB API response (usually YYYY-MM-DD).
None: If no data is found for the given geometry and GRB type combination.
"""
if border_distance > 0:
geometry = create_donut(geometry, border_distance)
bbox = get_bbox(geometry)
actual_url = (
GRB_FEATURE_URL + "/" + grb_type.upper() + "/items?"
"limit=" + str(limit) + "&crs=" + crs + "&bbox-crs=" + crs + "&bbox=" + bbox
)
update_dates = []
collection = get_collection(actual_url, limit)
if "features" not in collection:
return None
for c in collection["features"]:
if intersects(geometry, shape(c["geometry"])):
versiondate = datetime.strptime(
c["properties"][GRB_VERSION_DATE], DATE_FORMAT
).date()
update_dates.append(versiondate)
update_dates = sorted(update_dates, reverse=True)
if len(update_dates) > 0:
return update_dates[0]
return None
def get_collection_grb_actual(
geometry,
grb_type=GRBType.ADP,
partition=1000,
limit=DOWNLOAD_LIMIT,
crs=DEFAULT_CRS,
date_start=None,
date_end=None,
):
url = (
GRB_FEATURE_URL
+ "/"
+ grb_type.upper()
+ "/items?limit="
+ str(limit)
+ "&crs="
+ crs
)
if grb_type == GRBType.ADP:
name_reference_id = GRB_PARCEL_ID
elif grb_type == "gbg":
name_reference_id = GRB_BUILDING_ID
elif grb_type == GRBType.KNW:
name_reference_id = GRB_KNW_ID
else:
logging.warning(
f"type not implemented: {str(grb_type)} -->No reference-data loaded"
)
return {}, None
versiondate_filter = ""
versiondate_filter_start = ""
versiondate_filter_end = ""
if date_start is not None:
versiondate_filter_start = (
GRB_VERSION_DATE + ">" + date_start.strftime(DATE_FORMAT)
)
versiondate_filter = versiondate_filter_start
if date_end is not None:
versiondate_filter_end = GRB_VERSION_DATE + "<" + date_end.strftime(DATE_FORMAT)
versiondate_filter = versiondate_filter_end
if not (date_start is None and date_end is None):
versiondate_filter = versiondate_filter_start + " AND " + versiondate_filter_end
if versiondate_filter != "":
url = url + "&filter=" + versiondate_filter + "&filter-lang=cql-text"
collection = get_collection_by_partition(
url, geometry=geometry, partition=partition, limit=limit, crs=crs
)
return collection, name_reference_id
def get_collection_grb_fiscal_parcels(
geometry,
year=str(datetime.now().year),
partition=1000,
limit=DOWNLOAD_LIMIT,
crs=DEFAULT_CRS,
):
url = (
GRB_FISCAL_PARCELS_URL + "/Adpf" + year + "/items?"
"limit=" + str(limit) + "&crs=" + crs
)
return get_collection_by_partition(
url, geometry=geometry, partition=partition, limit=limit, crs=crs
)
def get_collection_grb_parcels_by_date(
geometry,
date,
partition=1000,
limit=DOWNLOAD_LIMIT,
crs=DEFAULT_CRS,
):
collection_year_after = get_collection_grb_fiscal_parcels(
year=str(date.year),
geometry=geometry,
partition=partition,
crs=crs,
)
# Filter on specific date: delete all features > specific_date
# TODO: experimental loader; unclear if we have to use "year-1 & year" OR if we have to use "year & year + 1"
collection_year_after_filtered = deepcopy(collection_year_after)
logging.debug(len(collection_year_after_filtered["features"]))
if (
"features" in collection_year_after_filtered
and len(collection_year_after_filtered["features"]) > 0
):
removed_features = []
for feature in collection_year_after_filtered["features"]:
versiondate = datetime.strptime(
feature["properties"][GRB_VERSION_DATE][:10], DATE_FORMAT
).date()
if versiondate > date:
removed_features.append(feature)
collection_year_after_filtered["features"].remove(feature)
logging.debug(len(collection_year_after_filtered["features"]))
# if no features are removed, return the full collection of year_after
if len(removed_features) == 0:
return collection_year_after
# if features are removed, search for the features in year before
collection_year_before = get_collection_grb_fiscal_parcels(
year=str(date.year - 1),
geometry=geometry,
partition=partition,
crs=crs,
)
kept_features = []
if "features" in collection_year_before and len(collection_year_before) > 0:
for feature in collection_year_before["features"]:
for rf in removed_features:
geom_feature = shape(feature["geometry"])
geom_removed_feature = shape(rf["geometry"])
if intersects(geom_feature, geom_removed_feature):
intersection = safe_intersection(geom_feature, geom_removed_feature)
if intersection.area > 1:
if feature not in kept_features:
kept_features.append(feature)
# search for intersection and check if it more than x%
# keep these features
# add them to
collection_specific_date = deepcopy(collection_year_after_filtered)
filtered_features = collection_year_after_filtered["features"]
specific_date_features = filtered_features + kept_features
logging.debug(len(specific_date_features))
collection_specific_date["features"] = specific_date_features
return collection_specific_date
def update_to_actual_grb(
featurecollection,
id_theme_fieldname,
base_formula_field=FORMULA_FIELD_NAME,
max_distance_for_actualisation=2,
feedback=None,
attributes=True,
):
"""
Function to update a thematic featurecollection to the most actual version of GRB.
Important to notice that the featurecollection needs a 'formula' for the base-alignment.
:param featurecollection: Thematic featurecollection
:param id_theme_fieldname: property-fieldname that states which property has to be used as unique ID
:param base_formula_field: Name of the property-field that holds the original/base formula of the geometry, that has to be compared with the actual formula.
:param max_distance_for_actualisation: Maximum relevant distance that is used to search and evaluate resulting geometries. All relevant distance between 0 and this max_distance are used to search, with a interval of 0.1m.
:param feedback: (default None): a QGIS feedback can be added to push all the logging to QGIS
:param attributes: (boolean, default=True): States of all original attributes has to be added to the result
:return: featurecollection
"""
logger = Logger(feedback)
# Load featurecollection into a shapely_dict:
dict_thematic = {}
dict_thematic_props = {}
last_version_date = datetime.now().date()
for feature in featurecollection["features"]:
id_theme = feature["properties"][id_theme_fieldname]
geom = shape(feature["geometry"])
# logger.feedback_debug("id theme: " + id_theme)
# logger.feedback_debug("geometry (wkt): " + geom.wkt)
dict_thematic[id_theme] = geom
dict_thematic_props[id_theme] = feature["properties"]
try:
base_formula_string = feature["properties"][base_formula_field]
dict_thematic_props[id_theme][BASE_FORMULA_FIELD_NAME] = base_formula_string
base_formula = json.loads(base_formula_string)
logger.feedback_debug("formula: " + str(base_formula))
except Exception:
raise Exception("Formula -attribute-field (json) cannot be loaded")
try:
logger.feedback_debug(str(dict_thematic_props[id_theme]))
if (
LAST_VERSION_DATE in base_formula
and base_formula[LAST_VERSION_DATE] is not None
and base_formula[LAST_VERSION_DATE] != ""
):
str_lvd = base_formula[LAST_VERSION_DATE]
lvd = datetime.strptime(str_lvd, DATE_FORMAT).date()
if lvd < last_version_date:
last_version_date = lvd
except:
raise Exception(f"Problem with {LAST_VERSION_DATE}")
datetime_start = last_version_date
datetime_end = datetime.now().date()
base_aligner_result = Aligner(feedback=feedback)
base_aligner_result.load_thematic_data(DictLoader(dict_thematic))
base_aligner_result.name_thematic_id = id_theme_fieldname
affected, unaffected = get_affected_by_grb_change(
dict_thematic=base_aligner_result.dict_thematic,
grb_type=GRBType.ADP,
date_start=datetime_start,
date_end=datetime_end,
one_by_one=False,
geometry_thematic_union=base_aligner_result.get_thematic_union(),
crs=base_aligner_result.CRS,
)
logger.feedback_info(
"Number of possible affected OE-thematic during timespan: " + str(len(affected))
)
if len(affected) == 0:
logger.feedback_info(
"No change detected in referencelayer during timespan. Script is finished"
)
return {}
logger.feedback_debug(str(datetime_start))
logger.feedback_debug(str(base_formula_field))
# Initiate a Aligner to reference thematic features to the actual borders
actual_aligner = Aligner(feedback=feedback, max_workers=None)
actual_aligner.load_thematic_data(
DictLoader(data_dict=dict_thematic, data_dict_properties=dict_thematic_props)
)
actual_aligner.load_reference_data(
GRBActualLoader(grb_type=GRBType.ADP, partition=1000, aligner=actual_aligner)
)
actual_aligner.relevant_distances = (
np.arange(0, max_distance_for_actualisation * 100, 10, dtype=int) / 100
)
dict_evaluated, prop_dictionary = actual_aligner.evaluate(
ids_to_evaluate=affected, base_formula_field=BASE_FORMULA_FIELD_NAME
)
return actual_aligner.get_results_as_geojson(
resulttype=AlignerResultType.EVALUATED_PREDICTIONS,
formula=True,
attributes=attributes,
)
class GRBActualLoader(GeoJsonLoader):
def __init__(self, grb_type: GRBType, aligner, partition: int = 1000):
super().__init__()
self.aligner = aligner
self.grb_type = grb_type
self.part = partition
self.data_dict_source["source"] = grb_type.value
self.versiondate_info = {"name": GRB_VERSION_DATE, "format": DATE_FORMAT}
def load_data(self):
if not self.aligner.dict_thematic:
raise ValueError("Thematic data not loaded")
geom_union = buffer_pos(
self.aligner.get_thematic_union(), GRB_MAX_REFERENCE_BUFFER
)
collection, id_property = get_collection_grb_actual(
grb_type=self.grb_type,
geometry=geom_union,
partition=self.part,
crs=self.aligner.CRS,
)
self.id_property = id_property
self.input = dict(collection)
self.data_dict_source[VERSION_DATE] = datetime.now().strftime(DATE_FORMAT)
self.aligner.logger.feedback_info(f"GRB downloaded: {self.grb_type}")
return super().load_data()
class GRBFiscalParcelLoader(GeoJsonLoader):
def __init__(self, year: str, aligner, partition=1000):
super().__init__(_input=None, id_property=GRB_PARCEL_ID)
self.aligner = aligner
self.year = year
self.part = partition
self.data_dict_source["source"] = "Adpf"
self.data_dict_source[VERSION_DATE] = datetime(int(year), 1, 1).strftime(
DATE_FORMAT
)
self.versiondate_info = {"name": GRB_VERSION_DATE, "format": datetime_format_TZ}
def load_data(self):
if not self.aligner.dict_thematic:
raise ValueError("Thematic data not loaded")
geom_union = buffer_pos(
self.aligner.get_thematic_union(), GRB_MAX_REFERENCE_BUFFER
)
collection = get_collection_grb_fiscal_parcels(
year=self.year,
geometry=geom_union,
partition=self.part,
crs=self.aligner.CRS,
)
self.input = dict(collection)
self.aligner.logger.feedback_info(f"Adpf downloaded for year: {self.year}")
return super().load_data()
class GRBSpecificDateParcelLoader(GeoJsonLoader):
def __init__(self, date, aligner, partition=1000):
logging.warning(
"Loader for GRB parcel-situation on specific date (experimental); Use it with care!!!"
)
try:
date = datetime.strptime(date, DATE_FORMAT).date()
if date.year >= datetime.now().year:
raise ValueError(
"The GRBSpecificDateParcelLoader can only be used for dates prior to the current year."
)
except Exception:
raise ValueError(
"No valid date, please provide a date in the format: " + DATE_FORMAT
)
super().__init__(_input=None, id_property=GRB_PARCEL_ID)
self.aligner = aligner
self.date = date
self.part = partition
self.data_dict_source["source"] = "Adp"
self.data_dict_source[VERSION_DATE] = date.strftime(DATE_FORMAT)
self.versiondate_info = {"name": GRB_VERSION_DATE, "format": datetime_format_TZ}
def load_data(self):
if not self.aligner.dict_thematic:
raise ValueError("Thematic data not loaded")
geom_union = buffer_pos(
self.aligner.get_thematic_union(), GRB_MAX_REFERENCE_BUFFER
)
collection = get_collection_grb_parcels_by_date(
date=self.date,
geometry=geom_union,
partition=self.part,
crs=self.aligner.CRS,
)
self.input = dict(collection)
self.aligner.logger.feedback_info(
f"Parcels downloaded for specific date: {self.date.strftime(DATE_FORMAT)}"
)
return super().load_data()