Fourth in a GIS 💖 AI mini-series. I was interested in possible links in Béarn, SW France where I am now, to Ley Lines I learned about when in England. Alfred Watkins reinterpreted in 1921-1925 possible alignments of monuments, tumuli & churches (BBC). Part I yielded only six recorded archaeological points: the Béarn exercise is a landscape case study rather than a statistical test of alignments; GIS can reveal shared relationships with relief, water, routes, and geology, but cannot demonstrate intentional alignment or validate ley-line claims. Here we turn to Cambridgeshire (East Anglia, where I was before) that has over 100 monuments & other data to expand this project.
Note: as in all my AI-assisted posts and AI-guided work, my text is in italics.
Quick take
There are no apparent strong lineaments or Ley Lines to map in central East Anglia from the evidence collected. Would casting a wider net from Béarn thru Cambridgeshire & Peterborough to country-wide help? That would tax the Google Colab (scripting) & Open Street Map or BGS GeoIndex Onshore (data sources) internet use & download; these are kept free under "fair use" for personal research. It took for example over 5 min. to retrieve OSM data for Cambridgeshire & Peterborough alone, and over-traffic errors started during EU afternoons... when US started to get online!
This is all under Creative Commons (CC BY-SA 4.O): the prompt thread beginning Part I is here for you to re-use; warning, there are not a few twists & turns, so drop a comment below if specifics interest you.
Cambridgeshire GIS Workflow
This workflow extracts OpenStreetMap (OSM) data for Cambridgeshire using a county-wide bounding box, then summarizes OSM line features within 1 km buffers around representative points.
1. Define the study area in QGIS
Input: ONS Cambridgeshire & Peterborough polygon layer.
Outputs:
cambs_sites— two representative pointscambs_buffers_1km— two 1 km buffers- A county-wide bounding box in latitude/longitude
The script creates one point inside each polygon, buffers each point by 1 km, and then calculates a bbox from the full extent of the ONS layer, not from the buffers.
The bbox is a rectangle covering the whole county area, not just a corridor.
2. Extract OSM data in Colab
Input: The county-wide bbox from QGIS.
Output: A GeoPackage containing OSM roads, paths, and waterways.
The Colab script downloads OSM line features within the county-wide bbox, clips them to that rectangle, separates roads from paths, and exports the result as a GeoPackage.
It creates three layers:
osm_roads— roads and similar highway linesosm_paths— footways, tracks, cycleways, etc.osm_water— rivers, streams, canals, drains, ditches
3. Download the GeoPackage
Output: A
.gpkg file saved to your computer.
After the extraction finishes, run the download cell. This saves the OSM GeoPackage from Colab to your local Downloads folder.
4. Load OSM data into QGIS
Input: The downloaded GeoPackage.
Output: OSM layers available for analysis.
Add the GeoPackage to QGIS and rename the roads layer to osm_roads if needed. Script 3 expects osm_roads, osm_paths, and osm_water.
5. Run the analysis in QGIS
Input: OSM layers plus
cambs_buffers_1km.Output: A summary table with OSM statistics per buffer.
Script 3 takes each 1 km buffer and calculates:
- How many OSM road, path, and waterway lines intersect it.
- Total length of those lines inside the buffer, in metres.
It produces two rows—one per buffer—with fields such as:
road_countroad_length_mpath_countpath_length_mwater_countwater_length_m
Key distinction
| Item | What it is | Where it comes from |
|---|---|---|
cambs_sites |
Two representative points | Derived from the two polygons in the ONS layer |
cambs_buffers_1km |
Two 1 km buffers | Created around the two representative points |
osm_roads |
OSM road lines | Downloaded from OSM within the county-wide bbox |
osm_paths |
OSM path lines | Downloaded from OSM within the county-wide bbox |
osm_water |
OSM waterway lines | Downloaded from OSM within the county-wide bbox |
cambs_summary_*.gpkg |
Buffer-level OSM summary | Created by script 3 |
Maps & table
![]() |
| click to enlarge, original here |
![]() |
| click to enlarge, original posted here |
Why the first attempt looked sparse
The intended study area covered Cambridgeshire and Peterborough, together forming the geographical footprint of the Cambridgeshire & Peterborough Combined Authority. Established in 2017, the Combined Authority brings together the area's councils rather than replacing them with a single merged council.
In the first attempt, however, the workflow used two representative centre points instead of the full administrative polygons. The initial bounding box was calculated from the combined extent of two 1 km buffers around those points. With the points arranged roughly north–south, this produced a narrow corridor rather than a rectangle covering the whole study area.
That explains why the OpenStreetMap extract initially looked sparse: the search covered only that restricted corridor. The apparent lack of features was a consequence of the search extent, not an OpenStreetMap retrieval error.
Rerunning the workflow with the full multi-polygon geometry produced the intended county-wide search extent. The two attempts therefore illustrate an important distinction: buffering representative points is not the same as using the areas those points represent.
Current setup
The updated script 1 still creates the two site points and buffers, but now calculates the bbox from the full extent of the ONS layer. This means the OSM extract covers a much larger area, while script 3 still summarizes OSM within the two 1 km buffers.
Monuments, routes and radon: mapping proximity without mixing the signals
- How close is each monument to a mapped road or path?
- Where do the recorded monuments cluster?
- What radon-potential category is mapped in the surrounding grid square?
The central distinction is proximity. Measuring how near a monument is to a road, path or mapped area does not, by itself, demonstrate directional alignment, historical intent or a causal relationship.
Building the monument dataset in QGIS
The starting point was a monument point layer, together with OpenStreetMap road and path layers. The working data were reprojected to British National Grid (EPSG:27700), allowing the nearest-feature calculations to use meters.
Two successive nearest-feature joins attached information about the nearest mapped road and nearest mapped path to each monument. The resulting point layer retained the monument attributes alongside the joined route attributes and distance fields.
These distances describe the shortest mapped separation between a monument and the relevant line feature. They are not walking distances, measures of public access, or evidence that the modern route existed when the monument was placed. Their interpretation also depends on the completeness and classification of the source mapping.
Add geologic faults & aeromagnetic data?
| click to enlarge, original here |
Adding radon context
The next step added the downloaded Radon Indicative Atlas as a polygon layer. Our area-of-interest extract contained 6,450 grid tiles, with values spanning classes 1 to 6 and most tiles in class 1.
Rather than calculating distance to these tiles, the workflow identified which tile intersected each monument point and attached its radon attributes. This is a containing-area comparison, not another nearest-line calculation.
The important source field is CLASS_MAX. It records the highest mapped radon-potential class occurring anywhere within a 1 km grid square. It does not necessarily describe the precise ground beneath the monument.
![]() |
| click to enlarge, original here |
These bands concern the estimated percentage of homes exceeding the indoor radon Action Level. They are not outdoor gas measurements, percentages of radioactive ground, or direct measurements at the monument locations. The atlas supplies regional context; it is not a fault map or an independent measurement of ground permeability.
What the new columns record
-
radon_TILE: the reference of the selected grid square. -
radon_CLASS_MAX: the maximum mapped radon-potential class within that square. -
radon_VERSION: the version recorded in the source dataset. -
radon_band: a readable percentage-band label. -
radon_matches: the number of tiles intersecting the point. -
radon_status: a flag showing whether the point had a single match, multiple matches, no match or another issue.
The updated script preserved one output row for every input alignment row. Where a point intersected more than one tile, it selected the highest valid class and flagged that choice. Unmatched points were retained rather than silently discarded.
A new GeoPackage stored the enriched point layer, leaving the earlier results untouched. Because it already contains the monument geometry and inherited alignment attributes, a further join back to the original monument layer is optional.
Keeping the maps interpretable
A heatmap weighted by radon class initially offered an eye-catching view, but it combined two different things: monument clustering and the numeric class assigned to each point. We therefore separated the displays.
1. Radon context with monument locations
The radon polygons are best shown with categorized colors using CLASS_MAX, following the six atlas bands. Ordinary monument points sit above them. This keeps the grid boundaries visible and makes clear that the radon information belongs to an area rather than to a measured point.
2. Monument density
A separate heatmap gives every monument equal weight. It answers the simpler question: where do the recorded monuments cluster?
The heatmap radius controls how far each point's influence spreads. It should be documented and kept consistent when comparing maps. Any duplicate monument records must also be checked before interpreting the density pattern.
3. Distance to roads and paths
For individual distances, graduated point symbols are clearer than a distance-weighted heatmap. Each monument can be colored or sized according to its own verified road-distance or path-distance field, with the relevant network visible underneath.
A distance-weighted heatmap would combine distance with clustering: several moderately distant monuments could create a stronger hotspot than one very remote monument. Graduated points avoid that ambiguity.
Here, “isolation” must always be qualified. Distance from a mapped road is not the same as distance from a mapped path, and neither measures isolation from other monuments.
Presenting the comparison
The radon-context map and the equally weighted monument-density map are placed side by side. Note:
- Radon classes describe containing grid squares, not measurements at monument locations.
- Monument density does not demonstrate association with radon.
![]() |
| click to enlarge, original here |
Checks before drawing conclusions
Matching input and output row counts is necessary, but it is not enough. Nearest-feature joins can return multiple equally close matches, so the number of unique monument identifiers should also be checked. The radon script preserves input rows; it does not automatically resolve duplicates inherited from earlier joins.
- Confirm the working coordinate reference system and distance units.
- Verify which distance field belongs to roads and which belongs to paths.
- Check unique monument identifiers as well as total row counts.
- Inspect unmatched, multiple-tile and invalid-class radon records.
- Remove temporary filters before calculating full-dataset totals.
- Use a persistent monument identifier for any later joins, rather than the GeoPackage row key or radon tile reference.
The next analytical question
Finding monuments in higher-class tiles is not enough to demonstrate an association. The relevant comparison is whether those tiles contain a larger share of monuments than their share of the study area would suggest.
The next step is therefore to summarize unique monuments by radon class, then compare their percentages with the percentage of the exact area of interest occupied by each class. Tiles crossing the study boundary should be clipped and compared by area rather than simply counted.
Even a difference between monument percentages and area percentages would remain exploratory. Settlement patterns, land use, survival, accessibility and recording effort could all influence the observed distribution. A formal test would need an appropriate comparison design.
For now, the result is a reusable spatial dataset and a clearer mapping method: show radon context as areas, monument clustering as density, and individual route proximity as distances. Keeping those signals separate makes the patterns easier to examine without turning an attractive map into an unsupported explanation.
Data and software
- QGIS for reprojection, nearest-feature joins, spatial matching, symbology and print layouts.
- OpenStreetMap source layers for the mapped monuments, roads and paths. © OpenStreetMap contributors.
- BGS–UKHSA Radon Indicative Atlas for the radon grid and class definitions. The source dataset's license and attribution requirements should accompany any published maps.
Reading the wider landscape
Visually, the radon-class pattern appears to pick out parts of the northern and central lowland landscape extending towards the Wash. This is a qualitative observation, not a calculated relationship. The radon atlas provides geological context rather than a substitute for a bedrock or elevation map.
The flood-risk and sea-level-rise maps below add separate perspectives on the landscape - see here for more on this topic - they are presented here for descriptive comparison; no additional proximity calculations or statistical tests have been applied to these layers.
![]() |
| click to enlarge, original here |
Each map should be read according to its own source, date, legend and assumptions. Flood-risk mapping and sea-level-rise scenarios should not be treated as interchangeable, or as reconstructions of the landscape when the monuments were established.
Concluding remark
As mentioned b re: faults & aeromagnetic data, the maps above most definitely show congruences at a regional scale. In the prompt thread, AI mentioned:
Call radon geological context, rather than a “bedrock substitute”. BGS derives radon potential from both bedrock and superficial geology, combined with indoor measurements. So the pattern you’re noticing may be useful, but it cannot independently establish bedrock boundaries or low-lying terrain.
So Radon definitely outlines the high ground running N-S to the E, as well as lesser one running SW-NE in the Center. This Space Shuttle view will give you an appreciation of the relief in East Anglia. Note that views are looking south and Cambridgeshire is always at top right.
![]() |
| click to enlarge, original here |
| L is map location, R is picture perspective above, The Wash is blue |






