How can I plot street density over a polar histogram in QGIS or another GIS program?

Street density over polar histogram

I want to create a graph that qualitatively captures some characteristics of a street layout (directionality, density). I have tried “polar histograms” in GIS, which are a step in the right direction, however come with some clear limitations. Because only orientation and not location of line segments is considered, they are always symmetrical, which doesn’t represent street layouts with clear preferred directions very well (see figure).
With the tool "spatial analyst" ("line density") I have already determined the density of the lines in a pixel. Does anyone know how to convert this raster file into a polar histogram in QGIS or with another programme?

[enter image description here](https://i.stack.imgur.com/2AEOS.png)

Is there any way to plot density of street location - as opposed to orientation - as function of angle from origin (center) over the same polar histogram seen in the figure?

thanks!

I found your question quite interesting and I think it makes for a nice pet project. I have some code at hand that I made for a similar problem that I adapted to your situation. Comments are provided in the code, so I think its mostly self-explanatory. When providing the png from your post as input, it generates the outputs shown below:

import numpy as np
import matplotlib.ticker as plticker
from matplotlib import pyplot as plt
# read image 
im = plt.imread('2AEOS.png')
# remove 4th column
im = im[:,:,0:3]
# take average over color channels
im = np.mean(im,axis=2)
# identify dimensions of image and center point
dimensions = np.array(np.shape(im))
origin = np.round(np.flipud(dimensions)/2)
# draw input image and center point
f = plt.figure(1)
plt.gray()
plt.imshow(im)
plt.plot(origin[0],origin[1],'rs') 
plt.show()
# calculate angle of every pixels vector from origin to x-axis
# np.arctan2 function takes coordinate args in order [y,x]
degrees = np.zeros(np.shape(im))
for y in range(dimensions[0]):
for x in range(dimensions[1]):
point           = np.array([x,y])
v_origin_point  = point-origin
degrees[y,x]    = np.degrees(np.arctan2(v_origin_point[1],v_origin_point[0]))
# choose bin number
bin_number = 40
lspace     = np.linspace(-180,180,bin_number+1)
vals       = np.zeros(bin_number)
for i in range(len(lspace)-1):
# identify lower and upper bound of bin
low  = lspace[i]
high = lspace[i+1]
# sum pixel intensities in "im" where "degrees" is within range low-high 
vals[i] = np.sum(im[(degrees>=low) & (degrees<high)])
# Calculate thetas and ratii for circular plot
N = bin_number
theta = np.linspace(np.pi, -np.pi, N, endpoint=False)
# SURFACE AREA OF BIN SCALES WITH VALUE IN vals
radii = np.sqrt(vals)
width = 2*np.pi*np.ones(N)/N
# Plot in matplotlib
g = plt.figure(2)
ax = plt.subplot(111, projection='polar')
bars = ax.bar(theta, radii, width=width, bottom=0.0)
# Choose spacing of xticks and remove yticks
loc = plticker.MultipleLocator(base=np.pi/4) 
ax.xaxis.set_major_locator(loc)
ax.set_yticks([])
plt.show()

Street layout with centre-point

Calculated polar histogram

The histogram is obtained by assigning each pixel the angle (in degrees) between its vector to the centre-point and x-axis like so

Degrees over domain as greyscale

and summing up the pixel intensities that it finds in an interval of angles (a pizza slice of the image if you will). Just note that the sum of pixel intensity relates to the surface area of the histogram bins, not their radius. You can adapt this by changing the relation to radii = vals.