# Calculating total shared boundary between polygons based on their attributes?

I tried to do this over the weekend but had no luck - even after searching this forum. I'm trying to calculate the frontage (length of the property bordering the street) of a number of properties. I've got a shapefile of Hobart, Australia which I've slimmed down for this question as below.

Private Parcel = Houses

Casement = Streets ``````#Download and unzip the file from https://drive.google.com/open?id=18yyrDSFz2jp0N3R0ehs80leVs9Qyv1wy
#libraries
library(tmap)
library(rgdal)
#Import
#Visualise the data
#I want to calculate the total length of each 'Private Parcel' polygon shared border with 'Casement'. I don't know how to do this.
``````

Can anyone tell me how I might be able to do this?

After reading @Spacedman’s comment I noticed that some properties are aggregating.

``````#Previous Code
library(rgeos)
library(rgdal)
houses = example[example\$CAD_TYPE1 == "Private Parcel",]
frontages = gIntersection(houses, streets,byid=c(TRUE,FALSE))
length(frontages)
dim(houses)
plot(houses[12,])
houses\$frontage = gLength(frontages, byid=TRUE)

#My test
plot(houses)
#As you can see, the frontages are being aggregated between some properties.
`````` ``````#Here's my 2nd attempt
library(tidyverse)
library(sf)
# intersection
# find out about the length of each line segment
int\$len = st_length(int)
#Plot
ggplot(data=house) + geom_sf(fill=house\$CID,alpha=0.5,col=NA) + geom_sf_label(aes(label=CID),size=2) + geom_sf(data=int,col=int\$CID,size=2)
`````` Two properties are still aggregating - 1271016 and 1271014.

• Welcome to GIS SE! We're a little different from other sites; this isn't a discussion forum but a Q&A site. Please check out our short tour to learn about our focussed Q&A format. What happens when you run the code that you have presented?
– PolyGeo
Sep 9, 2019 at 2:36
• Where did you look? All geometric operations on `sp` objects are handled by the `rgeos` package, and you can use `gIntersection(house, roads)` to get the intersection geometry - which for you should be the frontage. Loop over houses, and use `gLength` to get the line length. Some things can go wrong but try it and re-post when you get stuck. Sep 9, 2019 at 12:14

Split out the houses and the street features:

``````> houses = example[example\$CAD_TYPE1 == "Private Parcel",]
> streets = example[example\$CAD_TYPE1 == "Casement",]
``````

Compute the intersection of houses and streets, using `byid` to compute the intersection for each house for all the street lines:

``````> frontages = gIntersection(houses, streets,byid=c(TRUE,FALSE))
``````

That gives us the same number of frontage features as houses:

``````> length(frontages)
 27
> dim(houses)
 27 19
``````

Lets see if this makes sense - pick a feature at random, number 12:

``````> plot(houses[12,])
`````` Looks good. If all you want is the length, extract this using `gLength` (again using `byid` to get the length of each of the 27 features) and add to the `houses`:

``````> houses\$frontage = gLength(frontages, byid=TRUE)
> spplot(houses,"frontage")
`````` Looks believable - the houses on the corners have more frontage.

Not sure if this code will have problems if, for example, the roads and houses don't exactly overlie, or if some house features have no frontage, or if houses cross street lines resulting in point intersections. But this works on your example and some adjustment might be needed once you scale up to more data...