## 1. Introduction

#### 1.1. Context and Aim of the Work

The safety assessments for newly built roads can be carried out according to the guidelines of the design standards. In Italy, specifically, the Ministerial Decree dated 5.11.2001 [

1], referring back to the road classification of the Highway Code [

2], defines the constructive, technical and functional characteristics of each type of road; with regard to these features, it imposes prescriptions for the design and verification of the geometric elements and their composition, so that the road users’ movement takes place regularly and safely.

The verification of the design accuracy includes the road consistency examination. The road consistency is the characteristic of the infrastructure that makes it suitable to meet road user expectations and implies that the geometric successive elements are coordinated in order to produce harmonious and homogeneous driver performance and to not generate incorrect and, sometimes, unsafe driving behavior [

3,

4,

5].

The rules that guide road design should ensure that a driver is spontaneously induced to adopt a driving behavior that is congruent with the real characteristics of the road. It is the case of the “self-explaining roads” [

6,

7], which are defined as roads where drivers perceive and recognize the features of the infrastructure they are driving on and “instinctively” know how to behave [

8]. Since the driver’s behavior is directly related to the perception of the route he/she is following, a designer should make the road as readable as possible, so as to guide the user in the behavioral choices to be adopted, suitable to deal with even the most dangerous situations [

9,

10,

11].

Many studies have found that some design choices, violating the user’s expectations, can negatively affect his/her behavior [

12,

13,

14,

15], by inducing the disoriented driver to change his/her usual driving choices and leading him/her to have risky driving behavior.

The inadequacy of the technical characteristics of the road can be detected under several aspects, including, for example, unsuitable adhesion values of tire-asphalt pavement interface to guarantee the transfer of dynamic actions, geometric horizontal alignments design or conditions of the environment surrounding the infrastructure insufficient in ensuring a correct perception of the route, especially for occasional users. Design, construction, but above all maintenance problems are often the cause or contributing cause of high levels of accidents on many road networks.

Therefore, it is necessary to investigate in more detail what the role of the infrastructure may be in the occurrence of accidental events, in order to assess what type of countermeasures can be implemented to ensure higher safety standards.

In particular, in order to be conducted, evaluations related to the consistency of the technical characteristics of existing roads, require data that correctly represent the geometric characteristics of the infrastructures; such data are often not accessible, not updated or do not have the accuracy and precision required to perform the analysis to be carried out.

In light of the above, the identification of the plano-altimetric geometrical characteristics of the existing roads is still a topic of great interest in the field of civil and transportation engineering. Considering these problems, the objective of this research is looking for an automated, economical and fast procedure able to identify the elements of a road layout belonging to an existing and operating road network. In particular, the horizontal geometry—and, more specifically, the continuous trend of curvatures—is especially needed in order to perform the ex post application of some standard design rules and assess the coherence, consistency and homogeneity of the road infrastructure.

There are no shared and reliable methodologies that establish how this activity should be carried out, for the purpose of performing road safety analyses and for the maintenance of the infrastructure. Although much research has presented techniques and procedures suitable to obtaining useful data from surveys, no methods are available for the extraction of the reference axis (it means a continuous and coordinated line, composed by basic elements mathematically described, like in the design process of a new track) of an existing road section. Thus, the aim of this work and the proposed contribution to the state of the practice is to overcome the described lack.

#### 1.2. State of the Practice

In general, in order to overcome the lack of information about the existing road networks, a lot of energy has been spent in research aimed at obtaining the geometry of the compositional elements of the layout from different sources. Various attempts were made to obtain information regarding the horizontal and vertical alignment both from aerial and satellite georeferenced images, both from road maps of geographic information systems, and from data collected with positioning sensors on mobile systems and vehicles.

There are many road geometry relief techniques to date, which can be divided into static and dynamic ones, depending on the instrumentation used to obtain the data.

Among the static procedures, the following may be enumerated:

- -
Fully automated acquisition methods that extract horizontal curve data directly from road maps in GIS environment [

16,

17].

- -
Automatic identification of road geometry from digital vector data [

18].

Among the dynamic ones, on the other hand, it is worth mentioning:

Detection systems using mobile mapping system (M.M.S.) technology. These are widely appreciated both for their versatility and for their low operating costs. In general, M.M.S. are made up of vehicles equipped with different instruments, properly integrated with each other: a satellite navigation device, an Inertial Navigation System (INS) and an odometer [

19,

20,

21,

22]. In specific cases, the “Digital Highway Data Vehicle” (DHDV) has been used, an integrated system through which three Euler angles, the driving speed and the vehicle acceleration were measured to perform the survey in a rigorous way [

20,

23]. Another instrumented vehicle able to collect pavement condition and asset data and not only geometric information is the van automatic road analyzer (ARAN) [

24].

Global Navigation Satellite System (GNSS) with GPS receivers mounted on trains [

25,

26,

27] or vehicles travelling at almost constant speeds, instrumented with vertical gyroscopes and gyro compasses able to provide information on the vehicle positions (x, y, z coordinates) and orientations (angle of pitch, roll and yaw) [

28,

29,

30,

31,

32]. In some cases, vehicle positions have been recorded via smartphone [

33].

The described techniques are generally efficient and valid for acquiring a dataset but, once the necessary points for the definition of the road geometry have been collected, there comes the need to develop a rapid and economic procedure that allows to easily and precisely obtain the geometric definition of the plano-altimetric road layout that should be mathematically formalized through simple and continuous geometric elements.

The above-mentioned studies and the proposed techniques, in fact, are able to obtain points and to determine the coordinates of the track, but they do not manage to define a continuous curvature graph (or a continuous reference line), which is useful to conduct safety evaluation by means of the ex post application of the design verifications. Only in some cases, they can efficiently identify circular curves with their radii and straight lines, but they cannot represent the reference road axis in the traditional design way as a succession of straight lines, circular arcs and clothoids.

Over the years, multiple algorithms, able to identify in an automated way the circular arcs and the straight lines, have been built; the most successful method is that of the least squares fitting of the points of the road axis. The most widely used and easily implementable algorithm is the least squares circle fitting, able to estimate radius, angle at the center and arc length with high accuracy; as can be seen in

Figure 1, it can minimize the deviations between the computed curve and the original points [

34].

The least squares regression method is a technique aimed at determining an analytical function that approximates a set of values without necessarily passing through the points themselves. In fact, if the data come from experimental measurements and are, therefore, affected by measurement or instrument errors, or if the data quality is not very accurate (few significant figures), then it is appropriate to approximate to least squares rather than interpolate. In particular, the function found must be the one that minimizes the sum of the squares of the distances between the points and the plotted curve [

35].

As an implementation of the aforementioned algorithm, an iterative least squares circle fitting algorithm has been proposed, which iteratively searches for circumferences approximating the available data, changing the center and radius of the arc from time to time, until the smallest possible differences between the input data and the originated output are achieved [

28]. Another similar method, useful for searching the characteristic parameters of horizontal curves, generates from the calculation of each circle that can be drawn using all the triplets of non-aligned points, reducing iteratively the distance between the circle and the complete set of points, using a minimization method [

36].

In general, the algorithms developed through the least squares fitting method have tried to evaluate the reliability and accuracy of the results obtained, considering the errors due to the surveying GPS instrumentation [

32,

37,

38].

Other approaches defined by data collected from inertial measurement units fitted to road vehicles have also been used for radius measurement, including the kinematic method, the geometry method, the lateral acceleration method and the chord offset method [

20,

39]:

The kinematic method plans to use the longitudinal speed (

v) and centrifugal acceleration (

a) values, directly obtained from the dataset provided by the vehicle, according to the following formula:

The geometry method involves linear approximation of the vehicle’s trajectory between one measurement point and the next one—linearization made possible by the high data sampling frequency—and calculates the curve radius through the reciprocal of the curvature:

The lateral acceleration method uses the data of super-elevation, side friction factor and speed to calculate the corresponding radius:

According to the chord offset method, the chord length (

L) and the offset between the chord and the arc midpoint (

M) are measured through the vehicle’s trajectory; using these parameters, it is possible to calculate the radius of the curve:

In order to identify of the starting and ending points of the circular arcs; instead, the heading angle trends have been studied and their linear variation has been determined [

20,

39].

Finally, it is worth mentioning that several authors have used a completely different method from what was exposed so far to recreate the road geometry; this method involves the interpolation of raw data, acquired to identify the road layout, through polynomial curves and/or splines [

29,

30,

40,

41].

In contrast to some previous research, which are based on the B-spline approximation [

30,

42,

43], the proposed method utilizes a data fitting process. The least square regression is an optimization technique that allows to find a function, represented by an optimal curve, that is as close as possible to a set of data. Other similar procedures have been implemented in these years to automatically identify all horizontal geometric elements and produce recreated alignment geometry for existing roads and railways [

29,

44,

45].

## 2. Materials and Methods

As anticipated, the objective of the analysis entails in the definition of a methodology able to recognize the reference axis, which means the road layout composing elements as a succession of straight lines, circular arcs and clothoids. This is important and necessary for a more general purpose, which consists in the preventive safety evaluation offered by existing roads, through procedures based on the verification criteria of safety conditions referred to the road geometry.

It should be noted that, for example, on Italian roads in 2018, road accidents with injuries to persons were 172,553, with 3334 victims and 242,919 injured (see

Figure 2). Accidents and injuries are progressively decreasing on urban roads and freeways, while they are increasing on rural roads.

Given the high accident rate found on single-carriageway rural roads, which are more suitable for studying and understanding the behavior of users in relation to the characteristics of the road, it was decided to focus the research on single-carriageway highways outside urban perimeters.

Since the Italian state roads are spread throughout the national territory and about 90% of this network is managed by the company ANAS S.p.A., the choice of the road network to analyze has fallen on the roads managed by this company, shown in

Figure 3.

The company ANAS S.p.A. has shared with the authors its road graph, realized by means of high-performance surveys and composed of three-dimensional geographical coordinates. The network graph is, therefore, made up of georeferenced points, but it does not recognize or describe the geometric elements making up the road. Thus, it was necessary to design and develop an original procedure, subsequently implemented in a programming platform, aimed at identifying the different geometric elements that constitute the planimetric road layout.

Therefore, the research essentially consists of a new approach proposed with the aim to reconstruct the horizontal road alignment by using traditional design elements such as straight lines, circular arcs and transition curves.

The first difficulty met has been the processing of a huge amount of raw data, specifically the georeferenced points in space. The first research step has been looking into filtering and sorting the vertices of the road graph to obtain an ordered set of geographical coordinates to be connected.

Once the heading direction graph was constructed, it has been necessary to clean the data from repetition and inaccuracy and then smooth the resulted line through the application of a moving average filter, as the values were affected by measurement or instrument errors.

The Savitzky–Golay filter used is a digital low-pass filter that allows to make the data trend more continuous; it is applied to a series of data points in order to decrease the signal noise without deforming the signal. The subsets of consecutive data points are fitted using a low order polynomial with linear least square method, so the convolution of all the polynomials is then obtained [

46]. The data, having a set of

n {

x_{j},

y_{j}} points, where

j = 1, 2, …,

n, and

x is an independent variable, whereas

y is an observed value, can be represented with a set of m convolution coefficients

Ci; the effect of convolution can be expressed as a linear transformation:

Three inputs are necessary to execute a Savitzky–Golay filter: the signal (

x), the polynomial order (

k) and its frame length (

f); in

Figure 4, an example can be seen.

Consequently, a data mining effort has been made to write a code able to automatically approximate the heading direction signal and generate a new trend, through which it has been easier to search for the starting and ending points of the sections with zero curvature.

Subsequently, using the vertices falling within the starting and ending points of each constant curvature element, a least square regression has been implemented in a different code, in order to automatically find centers and radii of the circular arcs, analyzing from time to time the deviation standards.

A linear least square problem is [

47]:

where

${f}_{j}\left(x,r\right)$ is the “residual”,

The optimal values of the original variables

x,

r can be recovered from the formulae:

Once the geometrical parameters were known, two other programming codes have been written to construct in an automatic and rapid way the curvature graph and the design speed profile, according to the Italian road design guidelines.

This procedure will be described below and Highway n. 4 (S.S. n. 4) will be reported as a case study.

## 3. Results

The main obtained results refer to the application to the case study of the “Strada statale n. 4 Via Salaria”. This is an Italian state highway (HWY); it is a single carriageway highway for most of its route, presenting one lane per each travel direction, but it falls into the category of dual carriageway urban road in its first part (near to Rome city center). It, therefore, follows that this infrastructure is characterized by many various geometrical elements, changing gradually from urban highway to mountain road.

Once the graph of the road network has been acquired form the company ANAS, the 3D geographical information of the road of interest has been extrapolated, in terms of x, y and orthometric elevation coordinates. Thanks to the east and north coordinates in the WGS 84 UTM zone 32N Reference System, and the elevations in meters above sea level (m.a.s.l.), the planimetric and then altimetric layout of the road have been reconstructed.

Conventionally the direction of the distance has been taken from point A to point B, linking Rome to Ascoli Piceno, as it can be seen in

Figure 5:

The plan reconstruction, as a succession of straight and curved elements, has been carried out in stages through the implementation of a programming code.

Initially, a denser discretization of the axis was carried out, measuring the distances between two consecutive vertices and dividing the segments longer than 6 m. For each point the Cartesian coordinates, the elevation, the identification code, the azimuth of the segment linking the i-th point to the next one (see

Figure 6), the slope (m) and the y-intercept (q) of the straight line passing between the 2 above points, the partial length and, lastly, the distance have been tabulated, as shown in

Table 1.

Once the azimuth (or heading angle) of each segment that composes the graph has been calculated, its trend and the elevations have been represented in a diagram as a function of the distance, as shown in

Figure 7:

Subsequently, the Savitzky–Golay filter was applied to the heading angle trend for the purpose of smoothing the data, to increase their precision without distorting the signal tendency.

The results of this filtering process are shown in

Figure 8, for the highway section used as an example, from the distance 40 km to the distance 55 km:

Subsequently, again by means of the programming code, a further method to determine the starting and ending points of the straight lines has been developed, through the analysis of the heading angle. A straight line has constant heading all the way along it, so in the distance–heading direction diagram, it is possible to identify in each plateau the elements with zero curvature. On the contrary, the circular arcs are represented by linear variations of the heading angle; conventionally, positive slope of the heading angle trend corresponds to right-handed curves; the negative slope corresponds to left-handed curves. When there are linear variations with a positive slope immediately followed by variations with a negative slope, or vice versa, in which no section with a zero slope is interposed, the points of reverse curvature can be identified at the vertices of the cusps.

Circular arcs with a length shorter than 45 m were excluded from the analysis, because they can be too short to carry out a reliable analysis, and the length may depend on excessive signal variations, not effectively smoothed by the application of the filter described above.

Figure 9 shows, between the distance 40 km and 55 km, the identification of the starting points of the straight lines in red, the ending points in green and the points of reverse curvature in cyan.

Once the coordinates of the starting and ending points of the straight sections have been identified, the starting and ending points of the circular arcs have also been defined, dividing each linear variations of the heading angle into five equal parts. The search of the geometric parameters of the elements—i.e., the radii and the coordinates of the centers of the arcs, the length of the straight lines—has been carried out by analyzing only the vertices falling in the central 3/5 of the different sections. The remaining areas—between the straight lines and the arcs or between two successive arcs—have a curved trend and have not been investigated due to the possible presence of transition curves inside the path (or in any case they can belong to the gaps where the trajectory of the vehicle used for the survey has made a progressive steering to enter the circular curve).

Once the vertices of the graph, corresponding to each geometric element, were known, the following features for the straight lines have been defined, as shown in

Table 2: the length, the mean and the deviation of the azimuth, the starting and ending point in east and north coordinates.

The vertices of the graph that fall within the starting and ending points of a circular arcs were used to calculate the center and radius of each circumference, using a code that implements a least squares circle fitting algorithm.

Table 3 shows the east and north coordinates of the center and the value of the radii of the circular arcs, the deviations, the east and north coordinates of the starting and ending points of each element, the lengths of the arcs and the values of the angles at the center:

The research findings consist of a procedure composed by three sequential phases, as described below.

#### 3.1. Geometrizing Horizontal Alignment of an Existing Road Layout

As shown in

Figure 10, the Italian highway S.S. n. 4 has been completely geometrized by a succession of curves and straight lines; all the circumferences obtained by the previously described method, implemented in the original programming code, have been imported into the CAD (Computer-Aided Design) environment.

Figure 11 shows an example of a circumference computed by the circle fitting procedure implemented in the code, in which the points used for the least squares arc evaluation—and the center and radius coordinates definition obtained as a result—are represented in red (the curve obtained is in black):

To test the accuracy of the results, after calculating the characteristic parameters of the circular arcs and the straight lines through the least square fitting, a qualitative comparison between the project axis obtained and the road layout, obtained from GIS basemaps, has been done. Even if it is not possible to consider

Figure 12 as a quantitative validation, it can be noted that the two trends are overlapped in the entire distance considered.

#### 3.2. Curvature Graph

Once the calculation of the planimetric elements parameters has been completed, the curvature graph of the entire road layout has been subsequently constructed through another original code. Conventionally, the right-hand curves have been positioned above the reference (zero-curvature) line, assuming as positive the direction from point A to point B, while the left-hand curves have been positioned below. The transition zones, belonging neither to the straight lines nor to the circular arcs, have been treated as variable curvature elements, in order to connect elements having different curvatures.

Figure 13 shows the curvature graph for the segment taken as an example between the distance 40 km and 55 km:

#### 3.3. Design Speed Profile

Once the curvature graph has been defined, an additional code has been implemented to calculate the design speed profile. The calculation follows the indications of the Italian Ministerial Decree 5 November 2001 [

1] according to which:

The speed is constant all along the curve extension, and when the radius is less than a characteristic one (R_{2.5}), it is determined by some graphs in the guidelines.

On straight lines, on circular arcs with radius bigger than R_{2.5} and on clothoids, the design speed tends to the upper limit of the speed range; the acceleration spaces resulting from the exit from a circular curve and the deceleration spaces for the entrance to a curve are only limited to the elements considered.

The acceleration and deceleration values are 0.8 m/s^{2}.

It is assumed that the longitudinal slopes do not influence the design speed.

The transition distance

DT is the length in which the speed, according to the accepted theoretical model, passes from the value s

_{d1} to s

_{d2}, of two consecutive elements.

DT [m] is given by the following expression, where

$\mathsf{\Delta}s$ is the difference between s

_{d1} and s

_{d2} [km/h],

$\overline{s}$ is the average speed [km/h]:

The result of the calculation is shown in

Figure 14:

After performing the described three sequential phases, a complete representation of the geometric elements composing the planimetric road layout, as well as the theoretical design speed trend along the alignment, is made available.

This kind of representation is useful when performing the safety analyses and, first of all, the ex post application of the regulatory model for the design verification. In this way, the procedure allows us to identify the sections where in-depth assessments have to be focused: if there are road sections in which some theoretical design conditions are not met significantly; in fact, further investigation should be conducted, especially when these sections also show other critical conditions related to safety or functional indicators (accident rates, heavy traffic conditions, black spots and so on).

As a general remark, the proposed assessments can be useful to pursue strategical and managing actions aimed at planning any road upgrading activities and at improving road safety performances.

## 4. Conclusions

The research carried out was aimed at implementing an automated, fast and economical procedure to identify the horizontal geometry of existing two-lane rural roads, based on the data obtained from the road graph of the network. In particular, it has been necessary to search for methods and algorithms capable of identifying the reference axes, such as in the design process of new tracks, which have to be mathematically described as simple and as continuous geometric elements (circular curves, clothoids and straight lines).

Using semi-automated procedures, implemented in a programming platform, the following main steps have been developed:

Representation of the heading direction as a function of the distance, on the basis of the 3D spatial coordinates of the vertices of the road graph.

Application of a Savitzky–Golay filter to the heading angle trend, with the purpose of smoothing the data and in order to identify straight line starting and ending points and points of reverse curvature.

Analysis of the vertices falling within each element range and determination of the azimuth and length of the straight lines, of the radii and of the lengths of the circular arcs, by the application of a least square fitting procedure.

Identification of the transition zones between the constant curvature elements (treated as clothoids) to compose a continuous curvature diagram.

Once the curvature graph has been defined, calculation of the design speed profile by the implementation of an additional code.

The proposed method can be useful to analyze the technical characteristics of existing roads, especially in order to perform the ex post application of the regulatory models for the design verification. In particular, it allows us to perform road consistency analysis or to recognize other relevant effects in the road–user behavior interaction, such as the violation of:

All these analyses are intended to assess Minimum and maximum lengths of each element with constant curvature, in relation to the drivers’ correct perception of the road layout.

Correct succession of straight lines and circular arcs, or of two circular elements, in order to have gradual variations between their geometric characteristics.

Optical, dynamic and road users’ comfort criteria of the variable curvature elements.

All these analyses are intended to assess the safety conditions for existing and open to traffic roads. In fact, considering the high accident rate found on the selected kind of roads, it is important to better study and understand the influence of the characteristics of the road on the users’ driving ability. Thus, the presented procedures are suitable for researchers to evaluate the effects of geometric inconsistencies on driver behaviors, as well as for practitioners (road administrators, technicians and professionals) to more easily develop the assessments on existing networks.

The effective contribution of this research to the state of practice is due to the opportunity to extract the track geometry from surveys or databases. The curvature graph and the design speed profile can be rapidly obtained, especially if the raw data are properly stored; the entire analysis process (starting from the data filtering to the final diagrams) can be uploaded and performed in a few hours. The short time needed is suitable for an extensive and systematic application of the methodology.

In fact, it will be possible to implement a code capable of quickly verifying the compliance of existing road layouts with the current design standards, in order to identify the sections that most differ from an optimal geometric configuration.