QGIS – How to Anonymize Site Locations in Shapefile

coordinate systemqgisrSecurityshapefile

I'm putting together some training materials for a workshop on R for spatial data analysis in archaeology, and need to anonymize the true site locations (confidential information) while keeping the spatial integrity of the features within the site's bounding box.

My initial thought is to subtract the minimum X & Y values of the bounding box (all are in UTM coords) from all of the geometries — i.e., make the whole thing an arbitrary relative grid.

I have ESRI point and polygon shapefiles imported to R with sf, and some rasters loaded with terra but not sure how to go about doing the global spatial adjustment for either.

I also use QGIS, so some solution through there to pre-process the geometries before loading to R would work as well — but again, don't know how to do that.

Best Answer

There's a function in the maptools package called elide that is designed to shift and rotate spatial data. It can also flip, reflect, and rotate vector data, but it might be hard to do the equivalent operation on raster data in the same locations.

Disguising the real location of data is its raison d'etre:

Description:

     Methods for function ‘elide’ to translate and disguise coordinate
     placing in the real world.

At the moment it only works on sp class objects, but you can do a conversion step and use it.

Example, given an sf data frame with geometry:

> nc$geometry
Geometry set for 100 features 
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: -9386880 ymin: 4012991 xmax: -8399788 ymax: 4382079
Projected CRS: WGS 84 / Pseudo-Mercator
First 5 geometries:
MULTIPOLYGON (((-9069486 4332934, -9077066 4338...
MULTIPOLYGON (((-9043562 4351030, -9043652 4352...
MULTIPOLYGON (((-8956335 4334068, -8958566 4335...
MULTIPOLYGON (((-8461241 4344709, -8462173 4347...
MULTIPOLYGON (((-8595797 4333852, -8597683 4330...

I can create a shifted version:

> nc2 = st_as_sf(elide(as(nc, "Spatial"), shift=c(9069486,-4332934)))
> nc2$geometry
Geometry set for 100 features 
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: -317393.5 ymin: -319943.4 xmax: 669697.9 ymax: 49144.94
CRS:           NA
First 5 geometries:
MULTIPOLYGON (((-0.02079826 0.2879136, -7579.53...
MULTIPOLYGON (((25923.5 18095.96, 25834.31 2003...
MULTIPOLYGON (((113151.4 1133.67, 110920.2 2813...
MULTIPOLYGON (((608245.2 11774.57, 607312.6 142...
MULTIPOLYGON (((473689.5 918.4896, 471803.2 -27...

Note this has shifted it so that one particular vertex of the first polygon is at roughly (0,0). Pick your shift values however you wish.

This works with points, lines, and polygon sf objects.

For terra raster objects, you can use the shift function:

> r
class       : SpatRaster 
dimensions  : 90, 95, 1  (nrow, ncol, nlyr)
resolution  : 0.008333333, 0.008333333  (x, y)
extent      : 5.741667, 6.533333, 49.44167, 50.19167  (xmin, xmax, ymin, ymax)
coord. ref. : lon/lat WGS 84 (EPSG:4326) 
source      : elev.tif 
name        : elevation 
min value   :       141 
max value   :       547

Note change of extent when shifted:

> shift(r, 100,100)
class       : SpatRaster 
dimensions  : 90, 95, 1  (nrow, ncol, nlyr)
resolution  : 0.008333333, 0.008333333  (x, y)
extent      : 105.7417, 106.5333, 149.4417, 150.1917  (xmin, xmax, ymin, ymax)
coord. ref. : lon/lat WGS 84 (EPSG:4326) 
source      : elev.tif 
name        : elevation 
min value   :       141 
max value   :       547
Related Question