Mathematical formulation and preliminary testing of a spline approximation algorithm for the extraction of road alignments

Automation in Construction 47 (2014) 1–9

Contents lists available at ScienceDirect

Automation in Construction journal homepage: www.elsevier.com/locate/autcon

Mathematical formulation and preliminary testing of a spline approximation algorithm for the extraction of road alignments Laura Garach a,⁎, Juan de Oña a, Miguel Pasadas b a b

TRYSE Research Group, Department of Civil Engineering, University of Granada, Spain Department of Applied Mathematics, University of Granada, Spain

a r t i c l e

i n f o

Article history: Received 31 May 2012 Received in revised form 30 June 2014 Accepted 10 July 2014 Available online xxxx Keywords: Road alignments Highway design Segmentation Spline

a b s t r a c t This paper presents a methodology that allows the identification of the road alignments (curves, straights and clothoids) and their corresponding values of curvature based on a list of points with UTM coordinates obtained from field data. The procedure reconstructs the geometry of a road using a cubic spline that gives the road's singular points, geographically referenced, with the curvature values for each element. The methodology allows the user to select different parameters depending on the road type analyzed. It was applied in almost 1500 km of road with satisfactory results. The results of applying the methodology on a recently built road (whose alignments according to its project are known) are shown. A comparison of the project alignments with the alignments obtained according to the proposed method gives good results with small errors relative values, with maximum values of less than 4%. © 2014 Published by Elsevier B.V.

1. Introduction The horizontal geometry of a road is defined as a continuous axis formed by a succession of alignments or elements. Three types of alignments are normally used in road engineering: straight alignments in which the azimuth is constant and the curvature is null (infinite radius); circular curves in which the azimuth varies linearly along the trajectory and the curvature is constant; and transition curves, in which both the azimuth and the curvature vary along the trajectory [1]. The purpose of transition curves is to achieve a gradual transition in the change of direction from a straight to a curved alignment. The transition curve most frequently used for roads is the clothoid. The clothoid is a curve whose curvature changes linearly with its curve length and the product of a curve's radius by the distance to the point where the curvature is null is constant and equal to the clothoid squared parameter. Many applications used by transportation engineers require the road axis to be defined in terms of alignments, such as, for instance: maintenance and reconstruction works [2]; evaluation of driving conditions and the driver's mental workload [3]; analyses of road visibility [4]; road safety and operating speed; and design consistency [5,6]. A road's alignment is one of the geometric features that have the greatest impact on level of service and safety. Certain accidents frequency studies show that accidents on curved sections of roads are

⁎ Corresponding author at: ETSI Caminos, Canales y Puertos, c/Severo Ochoa, s/n, 18071 Granada, Spain. Tel.: +34 958 24 94 55. E-mail addresses: [email protected] (L. Garach), [email protected] (J. de Oña), [email protected] (M. Pasadas).

http://dx.doi.org/10.1016/j.autcon.2014.07.002 0926-5805/© 2014 Published by Elsevier B.V.

more frequent and severe than accidents on straight sections [7–10]. The determination of alignments is a preliminary step required by the Interactive Highway Safety Design Model [11] and the Highway Safety Manual [12] software (the former is used to analyze design consistency based on speed profiles and the latter provides a method to quantify changes in accident frequency as a function of cross-sectional features). The alignments of recently constructed roads are easy to obtain from design software. For existing roads, however, the data is not always available, or the format is not appropriate, it is not updated, or it is inaccurate. Significant progress has been made to obtain highway data with geospatial data capture technologies such as GPS (Global Positioning System), airborne and satellite-based remote sensing, map digitalization, and datalog vehicles. Most of these technologies provide a georeferenced set of points. A vehicle's trajectory or the centerline road axis can be reconstructed relatively easily using this data and this information suffices for many useful applications, such as navigation and route guidance systems. Previous research on centerline geometry extraction from geospatial data capture technologies, concentrates on laborious to fully automated systems that model road axis geometry in the form of a polyline or a line of continuous curvature [13,14]. Some researchers have focused on the use of various types of spline curves due to their flexibility, mathematical simplicity and computational ease [15–17]. A spline is a curve defined by a piecewise, explicit or parametric polynomial function. In approximation problems, polynomial spline functions are often used because they produce good approximation methods using low-grade polynomials and therefore prevent the undesired fluctuations caused by the polynomial approximation. The fact that low-grade polynomials usually perform better than those

2

L. Garach et al. / Automation in Construction 47 (2014) 1–9

of a high order, because the latter tend to induce high interpolation errors especially toward the data tail ends, is also applied in many engineering problems [18,19]. The most widely used spline functions for smoothness and order of approximation are odd degree splines as opposed to even degree splines and, among the former, cubic splines [17]. Class two cubic splines ensure the existence of a straight section that is tangent to every point on its diagram (class 1, continuous first derivative) and, on the other hand, the continuity of the curvature function all along its domain (class 2, continuous second derivative). Notwithstanding that splines can adequately describe centerline geometry, situations do exist where road axis geometry should be defined in terms of traditional design elements (i.e., straight lines, circular curves, and clothoids) [2,20]. Dividing the road into alignments depending on its curvature is a far more complicated task. In some works the identification of alignments is done manually, but this is not an accurate procedure and it is exceedingly tedious if used to analyze an entire road network. Today, the extraction of road alignments usually forms part of broader road inventory projects for which “land mobile mapping systems” are used for field work. Vehicle navigation for such systems relies on the use of integrated Global Positioning System (GPS), Inertial Measurement Unit (IMU) or Distance Measurement Instrument (DMI) through Kalman Filtering. Madeira et al. [21] and Chiang et al. [22], try to improve the Mobile Mapping Systems (MMS) technology used to obtain the road geometry, among other things. Wang et al. [23] recognize that MMS have yielded an enormous time saving in capturing road networks but considered that the manual extraction of the road information from the mobile mapping data is still a time-consuming task. Hence they develop a robust automatic road geometry extraction system developed by Absolute Mapping Solution Inc. (AMS). On the other hand, for railway or tunneling inventory projects, for which typically higher positioning accuracies are required, alternative vehicle navigation systems (such as those based on robotic total stations) may be used [24]. Several authors have proposed different sorts of algorithms to identify alignments based on data from geospatial data capture technologies [2,20,25–29]. Almost all these papers [2,20,25,26,28,29] identified two main phases. In the first stage, the georeferenced set of points is classified into straight lines, circular curves, and clothoids by comparing the angles of consecutive points [25,28,29] or by calculating and comparing the radius of consecutive points [2,20]. The second stage is a curve fitting stage in which the set of points is fitted to mathematical curves based on the analytical definition of the geometric elements (straight lines, circular curves, and clothoids): Jiménez et al. [20] use the minisum objective (minimization of the sum of the absolute value of distances between the set of points and the circumference) when fitting the alignments, whereas other authors [2,25,26,28,29] use a least squares regression approach. While most authors estimate the geometric parameters of alignments individually, contrary to existing methods that estimate the geometric parameters of highway and alignments individually, Tong et al. [29] discuss how to extract the entire set of parameters of combined alignments simultaneously using an integrated estimation method. In this paper we propose a method based on splines to identify road alignments (straight lines, circular curves, and clothoids), thus overcoming the main limitation that many authors have objected to in this technique [2,20]. Firstly, geospatial data capture technology is used to adjust a spline to the georeferenced set of points on the road. Then the spline is used to determine the curves, straight lines and clothoids. This method enables road alignments to be adjusted with or without clothoids. This proposal is based on the idea set out by Cafiso and Di Graziano [27]. However, they do not give a detailed description of or validate the methodology they propose, nor do they consider a solution for special cases, such as S curves.

This paper is divided into four sections. Following the Introduction, the proposed methodology is described in Section 2 and a case study in which the method is applied is presented in Section 3. The conclusions and a brief discussion for future research are presented in Section 4. 2. Methodology 2.1. The methodology A seven step methodology is proposed. The software used has been Mathematica which is a computer software system for doing mathematics [30]. 2.1.1. Step 1. Obtain the coordinates of a set of P points of the road The goal is to obtain the coordinates of a set of road points, normally uniformly spaced; using any type of geospatial data capture technology. Previously, the data must be pre-processed to filter errors in the data collection and eliminate any significantly outliers. Ben-Arieh et al. [16] develop a method for pre-processing data. Digital road maps can be made in several ways, by digitizing paper maps, taking aerial photographs or using a datalog vehicle. One of the most common techniques to obtain the road geometry is GPS positioning. In this paper let it be assumed that GPS data are used to obtain the road axis that is to be defined in terms of alignments. 2.1.2. Step 2. Adjustment of the variational cubic spline function s to the set of P points obtained Odd degree splines have graphical representations smoother than those of even degree. This fact together with the fact that low-grade polynomials usually perform better than those of high degree, led us to use cubical splines. Once we had chosen a cubic spline, class 2 (n − 1) was chosen to ensure continuous functions, and with continuous first and second derivatives (to ensure the continuity of the curvature function) as well. A class 2 variational cubic spline s that adjusts to the set of road P points whose coordinates were obtained is constructed. At this point in the process, a smoothing parameter intervenes, whose value determines a curve's degree of smoothing as opposed to the approximation of the road points. A spline is a piecewise polynomial function typically constructed using low order polynomial functions, jointed at breakpoints with certain smoothness conditions. The breakpoints are defined in this context as knots. Diminishing the distance between knots improves the approximation error until reaching a value after which the error remains constant. If n is the degree of the spline, in order to ensure the smoothness of the approximation, typically n − 1 continuity conditions should be fulfilled. Given a partition of [a, b] in m subintervals, Δm ¼ fa ¼ t 0 b t 1 b ⋯ b t m ¼ bg:

ð1Þ

S3(Δm) denotes the set of cubic spline functions of degree less than or equal to three and class C2(dim(S3(Δm)) = m + 3), i.e., every function s : [a, b] → R such that 2

i) s ∈ C ð½a; bÞ;

ii) s
½ti−1 ;t i  ∈ ℙð½t i−1 ; t i Þ;

∀i ¼ 1; …; m ;

where ℙ([ti − 1, ti]) is the space of all the restrictions of the polynomial functions of a degree less than or equal to three in the interval [ti − 1, ti] [31–33]. Let {B30, …, B3m + 2} the B-spline basis of S3(Δm).

L. Garach et al. / Automation in Construction 47 (2014) 1–9

Let {a0, …, an} ⊂ [a, b] and fP0 ; …; Pn g ⊂ ℝ2 be a point set given by their UTM coordinates. A search is made for a spline function s : ½a; b→ℝ2 such that sðai Þ ¼ Pi ; i ¼ 0; …; n; and sðt Þ ¼ ðs1 ðt Þ; s2 ðt ÞÞ with si(t) ∈ S3(Δm), i = 1, 2. Thus sðt Þ ¼

m þ2 X

3

αi Bi ðt Þ;

ð2Þ

i¼0

where α0, …, αm + 2 ∈ ℝ2 are the problem unknowns, which can be obtained as the solution of the minimization problem [34,35]. Find s ∈ ðS3 ðΔm ÞÞ2 such that 2

J ðsÞ ≤ J ðvÞ; ∀v ∈ ðS3 ðΔm ÞÞ ;

ð3Þ

3

2.1.3. Step 3. Compute the curvature function of the smoothing variational cubic spline s at a set of uniformly distributed points After compute the parametric smoothness variational cubic spline sðr Þ; r ∈ [0, 1], the arc parameter is constructed as the function l : [0, 1] → [0, L] defined by  r 0 lðr Þ ¼ ∫0 s ðt Þ dt; l ∈ ½0; 1;

ð6Þ

 1 0 where L ¼ ∫0 s ðt Þ dt. The function l(r), r ∈ [0, 1], is increasing; hence it exists the inverse function t : [0, L] → [0, 1] such that t(l(r)) = r, for all r ∈ [0, 1]. Thus, the smoothing variational cubic spline function can be considered as the function s : ½0; L→ℝ2 defined by

being J ðvÞ ¼

n X i¼0

sðdÞ ¼ sðt ðdÞÞ; d ∈ ½0; L:

D E2 ″ hvðaÞ−Pi i þ ε∫a v ðt Þ dt;

2

2

w21

b

w22 ;

ð4Þ

2

þ for any w ¼ ðw1 ; w2 Þ ∈ ℝ , and ε ∈ (0, + ∞), with hwi ¼ called the smoothness parameter. Thus, the result of this step is the computation of a parametric smoothness variational cubic spline of class C2 that approximates the given points of the road. Note that for this the following parameters are necessaries: • the parametric domain [a,b]; normally [a,b] = [0,1] is considered; • the number of subintervals of the partition of [a,b], this is the parameter m the dimension of the spline space S3(Δm) is m + 3; • the number of the approximation points, n, and the associated parameters {a0, …, an}; usually i D X

ai ¼

P j −P j−1

j¼1

n D X

P j −P j−1

E E ∈ ½0; 1;

i ¼ 0; …; n;

j¼1

is taken, • and finally, the smoothness parameter (ε).

ð7Þ

ð5Þ

Obviously the parameterized curves by sðt Þ, t ∈ [0, 1], and sðdÞ, d ∈ [0, L], are the same approximation curve of the given road. Now the curvature function k(d), d ∈ [0, 1] is considered, defined by D E 0 ″   s ðt ðdÞÞ  s ðt ðdÞÞ ds ; d ∈ ½0; L: ¼ kðdÞ ¼ dd hs0 ðt ðdÞÞi3

ð8Þ

Next, the curvature diagram of the smoothness variational cubic spline can be obtained, based on the arc parameter d (distance function) as shown in Fig. 1 (which represents a section of the A-4128 two-lane rural highway in the Province of Granada). 2.1.4. Step 4. Truncate the calculated curvature data so as to annul the values below a preset tolerance value The curvature diagram (see Fig. 1) shows a fluctuation in the curvature values owing to errors in the GPS data and the inevitable straying of the vehicle's trajectory from the road's axis. Therefore, a band of tolerance range taken from a constant width (± 1 / R) is defined, within which the curvature values are considered as nil, and anything inside this band is considered to be a straight section [27]. The value of the radius that defines the band depends on the type of road under study (a mountain road is not the same as a road on flat land).

-1

Curvature k (m )

Band of Tolerance 1/2000

Distance d (m) Fig. 1. Curvature diagram depending of the distance function and band of tolerance for filtering errors.

4

L. Garach et al. / Automation in Construction 47 (2014) 1–9

Having established the band of tolerance limited by k0 ¼ R10 , the truncated curvature function can be defined as the k∗(d), d ∈ [0, l], function given by



k ðdÞ ¼



Table 1 List of singular points with the correspondent curvature values.

kðdÞ if jkðdÞj N k0 ; 0; in otherwise;

whose graph is called a truncated curvature diagram (see Fig. 2). 2.1.5. Step 5. Determine the alignments The alignments are determined based on the truncated curvature function. To do so, under certain conditions, a trapezoid is adjusted to the relevant diagram of truncated curvatures at each non-null section of the diagram, which are determined by two consecutive null curvature points (see Fig. 2). The upper base of the trapezoid, characterized by a constant curvature value, enables the identification of curves and the sides of trapezoids in which the curvature varies from null for the straight sections to a non-null constant value. This enables the identification of clothoids. The software has been developed so that if the obtained clothoids length is less than a certain length that can be chosen by the user (e.g. 10 m), rectangles are adjusted instead of trapezoids (this is very useful for old roads where sometimes there are no clothoids). The corresponding trapezoid over the [di, di + 1], interval is uniquely determined by the following condition: the trapezoid's area is equal to the area of the region determined by the graph of the truncated curvature function and the abscissa axis between di and di + 1. Fig. 2 shows the application of the process of approximating the non-null sections of the diagram of truncated curvatures by trapezoids. Taking into consideration that the apexes of the polygonal determined by the trapezoids are singular points in the sense that they delimit the sections with different alignments (straight sections, curves and clothoids), the study ends with the list of coordinates for these singular points, the alignments determined by each pair of consecutive singular points and the values of the curvatures, as shown in Table 1.

UTM X

UTM Y

Distance

Curvature

Alignment

444998.4543 445073.8819 445125.3209 445251.2191 445290.1482 445356.3906 445390.7206 445427.7482 445468.5007 445507.0547 445524.2847 445582.8884 445591.0936

4148678.251 4148674.499 4148673.512 4148693.512 4148707.192 4148731.81 4148743.3 4148751.307 4148754.807 4148761.046 4148767.685 4148806.629 4148812.775

0 74.9879925 126.1075303 253.6505995 294.5919624 364.671598 400.6561106 438.4944427 479.2922693 518.9072725 537.7910773 607.8398094 617.9329431

0 0 0.002274778 0.002274778 0 0 −0.003816543 −0.003816543 0 0.007862313 0.007862313 0 0

Straight Clothoid Curve Clothoid Straight Clothoid Curve Clothoid Clothoid Curve Clothoid Straight Clothoid

2.1.6. Step 6. Correction of the truncated effect The truncating used in Step 4 causes the length of the straight alignments to increase and the length of the clothoids to diminish (see Fig. 3, which represents a section of the A-4128 two-lane rural highway in the Province of Granada). In step 4, the following alignments were identified: clothoid 1 of length CL1, straight of length SL, and clothoid 2 of length CL2. However, the actual length of the alignments should be CL1*, SL* and CL2*, respectively, where: CL1  ¼ CL1 þ L1

ð9:aÞ

CL2  ¼ CL2 þ L2

ð9:bÞ

SL ¼ SL−L1 −L2

ð9:cÞ

−1

ð9:dÞ

−1

ð9:eÞ

L1 ¼ ðR0 tagα1 Þ L2 ¼ ðR0 tagα2 Þ

where 1 / R0 is the value of the curvature in which the band of tolerance is established, tag is the tangent of the angle, and α1 and α2 are the

0.01 Truncated Curvature Trapezoids

Curvature k (m-1)

0.005

di

di+1

0 0

500

1000

1500

2000

2500

3000

-0.005

-0.01

Distance d (m)

Note: di,di+1 is the interval on which each trapezoid is adjusted Fig. 2. Approximation of the non-null sections of the diagram of truncated curvatures by trapezoids, to obtain alignments.

L. Garach et al. / Automation in Construction 47 (2014) 1–9

5

0.005 CL* 1

CL1

CL2

L2

L1

0.003

SL SL*

CURVE 1

Cuvature k (m-1)

CL* 2

Curvature Trapezoids

CLOTHOID 1

0.001

α1

1/R0

-1/R0

-0.001

BAND OF TOLERANCE

CLOTHOID 2

α2

-0.003 CURVE 2

Distance d (m)

-0.005 0

100

200

300

400

500

600

700

800

Note: SL*:Real Straight Length; CLi*:Real Clothoid i Length; SL: Straight Length according to Truncated; CLi:Clothoid i Length according to Truncated; Li: Distance belonging to the clothoid i which is added to the straight section by the truncating process Fig. 3. Correction of the truncated effect.

angles that form clothoids 1 and 2, respectively, with the horizontal line (see Fig. 3). Owing to the truncating lengths L1 and L2 are treated in step 4 as though they were part of the straight, when in fact they are part of the clothoids of each preceding and following curves. In step 6, this is corrected so the length of the final straight considered will be SL*.

This correction also permits a solution to the case of S curves with no intermediate tangent. The result of the truncating in step 4 is that when S curves exist, a straight line that is SL (= L1 + L2) long is always identified, when in fact it does not exist (see Fig. 4, which represents a section of the A-4128 two-lane rural highway in the Province of Granada). When the truncating is corrected, L1 and L2 are considered

0.0050 CL2*

CL1*

0.0040 CL1

L1 L2

CL2

SL

SL*=0

0.0030 CURVE 1

0.0020

CLOTHOID 1

Cuvature k (m-1)

0.0010

α1

1/R0

0.0000 -1/R0

-0.0010

α2

BAND OF TOLERANCE

CLOTHOID 2

-0.0020 -0.0030

CURVE 2

-0.0040 Curvature Trapezoids

-0.0050 -0.0060 0

50

100

150

200

250

300

350

400

450

Distance (m) Note: SL*: Real Straight Length; CLi*:Real Clothoid i Length; SL: Straight Length according to Truncated; CLi: Clothoid i Length according to Truncated; Li: Distance belonging to the clothoid i which is added to the straight section by the truncating process Fig. 4. Correction of the truncated effect for S curves.

6

L. Garach et al. / Automation in Construction 47 (2014) 1–9

Table 2 Resulting parameters from the application of the program in 1432 km of roads. Curvature radius (m)

Smoothing parameter

Band of tolerance (1 / R) (m−1)

[20–100] [100–300] N300

10−9 10−8 10−7

1/1500 1/2000 1/3500

as parts of clothoids, and the “false straight line” of length SL that was identified between the two S curves disappears, in such a way that the length of the real straight line between the two curves is nil (SL* = 0). 2.1.7. Step 7. Recalculation of alignments affected by correcting the truncation Having corrected the error introduced by the effect of the truncation, the alignments are recalculated in the manner described in Step 5, taking the corrected lengths into account (real length of straight lines SL* and clothoids CL1* and CL2*).

may be adequate for road safety, since curves with a radius of more than 1000 m normally do not pose a safety hazard. Nonetheless, he specifies that too low a cut-off value could omit an important part of many curves and, to the contrary, values that are too high could create ‘false curve’ records. On the roads studied in Granada, for example, in sinuous sections with very small radius curves (less than 100 m), the straight sections were considered to be any curve with radius larger than 1500 m. However, in sections where the majority of the radii of the curves were relatively large (more than 300 m), the limit between the curve and the straight section was considered to be 3500 m. Selecting the adequate parameters, the application of the program on around 1500 km of roads the Province of Granada (Spain), the results were satisfactory. Table 2 shows the parameters used based on the radius of section's curvature. 3. Case of study 3.1. Data

2.2. Discussion on the parameters used The proposed method permits the identification of alignments in existing road. Software has been developed for the proposed methodology, in which the user can select different parameters depending on the type of road analyzed. Some of these parameters are the spline's smoothing parameter (ε) or the limit between curves and straight sections. Being able to select ε is of paramount importance. For example, in this case, it was found that when roads with curves all along their layout with not too small radii (more than 300 m) were analyzed, an ε = 10−7 in the spline gave an almost perfect adjustment between the road's real geometry and the spline geometry. However, this ε in very sinuous roads, with curves that have very small radii, gave inadmissible separations between the two layout values. In very winding roads, very precise results were again obtained, by using an ε = 10−9 instead of 10−7, in which the spline sacrifices smoothing at the cost of interpolation. In addition, users can select the band of tolerance (or limit between the curves and straights). The limit is established at different values, depending on the authors. In Spain, according to the National Geometric Design Standards [36], the limit between curves and straight sections is a 3500 m radius for conventional roads and 5000 m for multilane highways and freeways. In New Zealand, for instance, any road elements with a radius of more than 3000 m are considered as straight sections [37]. However, Koorey [37] suggests that a radius of 1000 m

The procedure was applied to obtain the alignments of 1432 km of roads in the Province of Granada (Spain). As an example, Fig. 5 shows the alignments in a section of the A-4134 (0.62 km) road and in a section of the A-4005 road (1.08 km), respectively. In Section 3.4 we have verified the potential of the proposed method in existing roads. As a case study we chose a recently built road (the A-348 road) because is easy to obtain the project plans for new roads. The Andalusian regional government gave us the plans, thus enabling us to compare the alignments (straight segments, curves and clothoids) defined in the project with the alignments identified according to the proposed method. The length of the road section studied is 5.1 km and it is located in a mountainous area. 3.2. Application of the method to the case of study Firstly, the UTM coordinates for the A-348 road's axis, spaced at 10 m from each other, were obtained from the geometric design lists. Next, a parametric adjustment variational cubic spline of class 2 is built (which will allow the horizontal geometry of a road to be obtained), based on the points obtained and the following parametric values from the problem: • Domain: [a, b] = [0, 1]. • Distance between knots: 40 m.

Fig. 5. Alignments in two existing road sections.

7

Curvature k (m-1)

L. Garach et al. / Automation in Construction 47 (2014) 1–9

Distance d (m) Fig. 6. A-348 road: curvature diagram of the approximation spline function in points uniformly distributed every 5 m.

• Number of polynomial sections: m ¼ Length of40the road : • Adjustment parameter: ε = 10−7. Next, the curvature function in a set of spline points evenly spaced at 5 m from each other is calculated and the curvature diagram shown in Fig. 6 is obtained. Following the indications of step 4 of the general methodology, a 1 m−1 (R0 = band of tolerance limited by the curvature k0 ¼ 3500 3500 m) is set. Therefore, it is considered that all the elements with a 1 m−1 are going to be straight. curvature of less than k0 ¼ 3500 Fig. 7 shows the truncated curvature diagram. Finally, the alignments are obtained. To do so, a trapezoid is adjusted for each spline section between two null curvature points, with the condition that the trapezoid's area is the same as the spline section's area (see Fig. 7). A software has been developed to implement the method for determining the alignments described in this paper and the program's output

shows the list of singular points coordinates (UTM X and UTM Y coordinates), the alignments determined (curve, straight or clothoid) and the corresponding curvature values, as well as the distance to the road's origin, enabling the road alignments to be drawn (see Fig. 8). 3.3. Error analysis Table 3 compares the alignments as per the project design and the alignments obtained with the proposed method. It can be observed that the successions of alignments in the design and in the model proposed are very similar. For example, according to the A-308 project, the sequence of alignments from the beginning of the road to the end of the first kilometer (641.95 cumulative distance) is Curve, Clothoid, Straight and Clothoid. If the alignments' sequence obtained by the method proposed in this paper is analyzed, it can be observed that the type of alignments obtained is the same: Curve, Clothoid, Straight and Clothoid. Moreover,

Truncated Curvature

0.004

Trapezoids

Curvature (m-1)

0.002

0.000

0

1000

2000

3000

4000

5000

6000

-0.002

-0.004

-0.006

Distance (m)

Fig. 7. A-348 road: truncated curvature diagram for a band of tolerance limited by R0 = 3500 m and approximation of the non-null sections of the diagram of truncated curvatures by trapezoids.

8

L. Garach et al. / Automation in Construction 47 (2014) 1–9

¯

Pk 5+000

Two Lane Rural A-348 Pk 0+000

the purpose of knowing the alignments is to try to relate their characteristics (such as the radius of the curve) with the occurrence of traffic crashes. 3.4. Testing the method on existing road

0

ALIGNMENTS

0.5

1

1.5 Kms

The alignments for existing roads are usually not available (either there is no access to the project plans or they do not exist). Therefore, it is not easy to check if the alignments obtained with the proposed method are correct. In two existing road segments (A-4131 and A-4005), the measures of several radii have been obtained by the proposed method and have been compared with the measures obtained by AutoCAD (a drawing program). Fig. 5 shows the value of these radii obtained by AutoCAD. Table 4 shows the comparison of the measures obtained en several curves of two existing roads by the proposed method and by AutoCAD. The relative error of the circular curves' radii is studied considering the error as the quotient between the difference in our method's value and the AutoCAD's value divided by the AutoCAD's value. The relative errors values obtained (with maximum values of less than 10%) are small.

CLOTHOID CURVE !

STRAIGHT

!

!

!

!

!

!

!

!

!

!

0

50

100

200

300 Ms

Fig. 8. Alignments on the A-348 road layout.

4. Conclusions if the relative error of the circular curves' radii is studied, considering the error as the quotient between the difference in our method's value and the value taken as exact (the project's value) divided by the project's value, very small relative errors values (with maximum values of less than 4%) are obtained. For example, the maximum difference between the project radius (260 m) and the proposed method radius (251.82) is obtained in the curve whose cumulative distance is 4036.16. This error is totally acceptable for certain areas such as road safety, where

A method has been developed that adjusts splines to field data and uses them to identify and reference road alignments (curves, straight sections and clothoids). The method enables their corresponding curvature values to be calculated quickly and easily for an entire road system involving thousands of kilometers. Most existing methods for identifying alignments are based on the classification of points according to their curvature and the subsequent

Table 3 Comparison of alignments according to the project and according to the proposed model. Cumulative distance

0.00 234.13 352.01 641.95 745.94 967.84 1091.83 1146.83 1238.83 1454.71 1566.71 1631.71 1691.40 2121.88 2206.57 2294.94 2548.06 2656.43 2720.23 3522.42 3601.21 3684.50 3857.84 3941.12 3951.13 4036.16 4125.96 4223.99 4303.26 4518.56 4620.82 4845.80 4928.82 5122.55

Project

Proposed method

UTM X

UTM Y

Alignment

494965.01 495153.38 495349.98 495569.17 495717.60 495895.53 496016.92 496032.59 496133.64 496319.03 496443.68 496482.26 496573.59 496949.86 497018.60 497087.06 497330.36 497417.94 497507.05 498193.58 498245.22 498295.58 498418.11 498419.25 498494.95 498569.15 498635.22 498670.98 498704.70 498888.64 498970.32 499212.95 499289.30 499443.41

4089370.20 4089343.65 4089355.11 4089376.82 4089397.79 4089469.76 4089557.97 4089570.39 4089643.82 4089710.81 4089718.91 4089719.81 4089723.93 4089875.31 4089935.59 4089993.15 4090054.56 4090036.39 4090015.95 4090284.22 4090359.67 4090432.04 4090543.42 4090544.15 4090586.66 4090631.00 4090704.94 4090783.63 4090858.04 4090978.71 4090980.02 4090970.67 4090962.81 4090775.73

Curve Clothoid Straight Clothoid Curve Clothoid Straight Clothoid Curve Clothoid Straight Clothoid Curve Clothoid Clothoid Curve Clothoid Clothoid Curve Clothoid Clothoid Curve Clothoid Straight Clothoid Curve Clothoid Clothoid Curve Clothoid Straight Clothoid Curve Clothoid

Rm = method radius; Rp = project radius. ⁎ Error = relative error radius.

Radius 811.94

600.00

−500.00

700.00

−350.00

700.00

−500.00

260.00

−250.00

−200.00

UTM X

UTM Y

Alignment

494958.31 495189.91 495307.55 495596.12 495699.09 495904.22 496004.45 496047.55 496121.90 496322.59 496434.29 496499.29 496558.92 496953.20 497014.26 497084.27 497323.35 497429.45 497491.58 498199.63 498246.63 498291.18 498417.90

4089372.04 4089343.60 4089350.98 4089379.50 4089393.90 4089475.17 4089548.08 4089582.26 4089636.43 4089711.27 4089718.73 4089720.18 4089722.81 4089878.07 4089931.96 4089991.03 4090055.22 4090033.54 4090019.01 4090292.49 4090361.87 4090426.24 4090543.19

Curve Clothoid Straight Clothoid Curve Clothoid Straight Clothoid Curve Clothoid Straight Clothoid Curve Clothoid Clothoid Curve Clothoid Clothoid Curve Clothoid Clothoid Curve Clothoid

498494.83 498571.72 498631.60 498671.65 498706.09 498878.64 498980.70 499205.57 499288.10 499429.12

4090586.51 4090633.18 4090699.32 4090785.39 4090860.01 4090977.08 4090979.78 4090971.17 4090962.84 4090841.79

Clothoid Curve Clothoid Clothoid Curve Clothoid Straight Clothoid Curve Clothoid

Radius

Error⁎ (Rm − Rp) / Rp

805.23

0.83%

595.28

0.79%

−495.49

0.90%

697.05

0.42%

−346.22

1.08%

699.50

0.07%

−492.07

1.59%

251.82

3.15%

−246.61

1.36%

−196.02

1.99%

L. Garach et al. / Automation in Construction 47 (2014) 1–9 Table 4 Comparison of radius according to AutoCAD and according to the proposed model. Road

Proposed method radius (Rm)

AutoCAD radius (Ra)

Error⁎ (Rm − Ra) / Ra

A-4134 Curve 1 A-4134 Curve 2 A-4005 Curve 1 A-4005 Curve 2

31.05 42.48 93.72 162.28

30.22 38.91 87.5 169.7

2.75% 9.18% 7.11% 4.37%

⁎ Error = relative error radius.

adjustment of geometric curves to the experimental data [20,25,26,29]. Although many such works demonstrate the mathematical and computational simplicity of splines and their high degree of flexibility [2,20], using them is discarded owing to the impossibility of identifying curves, straight sections and clothoids. The method proposed in this paper profits from the benefits of splines and also allows alignments to be identified, thereby overcoming the above-mentioned limitation. The computational simplicity of splines has enabled the methodology to be implemented on around 1500 km of existing two-lane rural roads in Spain. In addition, the software developed with the proposed method allows users to select the right parameters to adapt to different types of roads, such as the spline's smoothing parameter (ε) and the band of tolerance (1 / R). The result of the coordinates obtained can be the input of highly important road safety software, including the IHSMD and the HSM, which need the road alignments as baseline data. The proposed methodology provides a solution to most road situations, such as curves with or without clothoids (in old roads) and S curves. However, the method presented in this paper does not resolve the problem of ovoid curves. The authors will attempt to expand the algorithm to solve this particular case in future works. Another line of work planned for the future is to include the verification of specified standards of the geometry for each road under study (i.e. minimum lengths for straight sections, curves and clothoids) in the software. This will enable the software to detect lack of compliance with regulations in old roads built before the standards entered into force and to verify the compliance of new roads built after the regulation. The program will mark any sections that do not comply with regulations so they can be studied by the analysts as potential points of inconsistencies in the geometry or, in other words, they can be studied as points potentially connected to accidents, from the viewpoint of road safety.

Acknowledgments The researchers wish to thank the Andalusian regional government. Without the data they provided, this research would not have been possible. The authors appreciate the reviewers' comments and effort in order to improve the paper.

References [1] C. Kraemer, J.M. Pardillo, S. Rocci, M.G. Romana, V. Sánchez-Blanco, M. Del Val, Ingeniería de Carreteras, vol. IMac Graw Hill, Madrid, 2003. [2] V. Gikas, J. Stratakos, A novel geodetic engineering method for accurate and automated road/railway centerline geometry extraction based on the bearing diagram and fractal behavior, IEEE Trans. Intell. Transp. Syst. 13 (1) (2012) 115–126. [3] O. Pellegrino, An analysis of the effect of roadway design on driver's workload, Balt. J. Road Bridg. Eng. 4 (2) (2009) 45–53 (+Ia + IIa + IIIa). [4] M. Castro, L. Iglesias, J.A. Sánchez, L. Ambrosio, Sight distance analysis of highways using GIS tools, Transp. Res. C Emerg. Technol. 19 (6) (2011) 997–1005.

9

[5] S. Cafiso, A. Di Graziano, G. Di Silvestro, G. La Cava, Safety performance indicators for local rural roads: a comprehensive procedure from low-cost data survey to accident prediction model, Transportation Research Board 87th Annual Meeting. Washington, D.C., 2008. [6] A.M.P. Zuriaga, A.G. García, F.J.C. Torregrosa, P. D'Attoma, Modeling operating speed and deceleration on two-lane rural roads with Global Positioning System data, Transp. Res. Rec. 2171 (2010) 11–20. [7] G.M. Gibreel, S.M. Easa, I.A. El-Dimeery, Prediction of operating speed on threedimensional highway alignments, J. Transp. Eng. 127 (1) (2001) 21–30. [8] B.-J. Park, K. Fitzpatrick, D. Lord, Evaluating the effects of freeway design elements on safety, Transp. Res. Rec. 2195 (2010) 58–69. [9] G. Rey, D. Clair, M. Fogli, F. Bernardin, Probabilistic assessment of roadway departure risk in a curve, Veh. Syst. Dyn. 49 (10) (2011) 1649–1671. [10] Y. Wang, E. Wenjuan, D. Tian, G. Lu, G. Yu, Highway curve speed precaution system based on VII. ICTIS 2011: multimodal approach to sustained transportation system development — information, technology, implementation, Proceedings of the 1st Int. Conf. on Transportation Information and Safety, 2011, pp. 1189–1196. [11] Federal Highway Administration (FHWA), IHSDM Design Consistency Module. Engineer's Manual, FHWA, Washington, 2004. [12] Federal Highway Administration (FHWA), Highway Safety Manual (HSM), 2011. [13] A. Drakopoulos, E. Örnek, Use of vehicle-collected data to calculate existing roadway geometry, J. Transp. Eng. 126 (2) (2000) 154–160. [14] L.B. Cremean, R.M. Murray, Model-based estimation of off highway road geometry using single-axis LADAR and inertial sensing, Proceedings of IEEE International Conference on Robotics and Automation. Orlando, Florida, 2006, pp. 1661–1666, (Article number 1641945). [15] Gruen, Li, Semiautomatic linear feature extraction by dynamic programming and LSB-Snakes, Photogramm. Eng. Remote Sens. 63 (8) (1997) 985–995. [16] D. Ben-Arieh, S. Chang, M. Rys, G. Zhang, Geometric modeling of highways using Global Positioning System data and B-spline approximation, J. Transp. Eng. 130 (5) (2004) 632–636. [17] M. Castro, L. Iglesias, R. Rodríguez-Solano, J.A. Sánchez, Geometric modelling of highways using Global Positioning System (GPS) data and spline approximation, Transp. Res. C Emerg. Technol. 14 (4) (2006) 233–243. [18] N. Sneeuw, Global spherical harmonic analysis by least-squares and numerical quadrature methods in historical perspective, Geophys. J. Int. 118 (1994) 707–716. [19] V. Gikas, P.A. Cross, A. Akuamoa, A rigorous and integrated approach to hydrophone and source positioning during multi-streamer offshore seismic exploration, Hydrogr. J. 77 (1995) 11–24. [20] F. Jiménez, F. Aparicio, G. Estrada, Measurement uncertainty determination and curve-fitting algorithms for development of accurate digital maps for advanced driver assistance systems, Transp. Res. C 17 (3) (2009) 225–239. [21] S. Madeira, J. Gonçalves, L. Bastos, Sensor integration in a low cost land Mobile Mapping System, Sensors 12 (3) (2012) 2935–2953. [22] K. Chiang, K. Bang, Y. Tseng, Y. Li, Direct sensor orientation of a land-based Mobile Mapping System, Sensors 11 (7) (2011) 7243–7261. [23] C. Wang, T. Hassan, N. El-Sheimy, M. Lavigne, Automatic road vector extraction for Mobile Mapping Systems, Int. Arch. Photogramm. Remote. Sens. Spat. Inf. Sci. XXXVII (2008) 516–522 (Part B3b, Beijing, China). [24] V. Gikas, S. Daskalakis, Determining rail track axis geometry using satellite and terrestrial geodetic data, Surv. Rev. 40 (2008) 392–405. [25] M. Imran, Y. Hassan, D. Patterson, GPS–GIS-based procedure for tracking vehicle path on horizontal alignments, Comput. Aided Civ. Infrastruct. Eng. 21 (5) (2006) 383–394. [26] S. Choi, J. Sung, Data generalization algorithm for extraction of road horizontal alignment design using the GPS/INS data in Lectures in Computer Science, Advances in Hybrid Information Technology, vol. 4413, Springer, Berlin/Heidelberg, 2007. 51–62. [27] S. Cafiso, A. Di Graziano, Automated in-vehicle data collection and treatment for existing roadway alignment, 4th International Gulf Conference on Roads Compendium of Papers, CD-ROM. Doha, Qatar, 2008. [28] A. Karamanou, K. Papazisi, D. Paradissis, B. Psarianos, Precise estimation of road horizontal and vertical geometric features using mobile mapping techniques, 6th International Symposium on Mobile Mapping Technology, Sao Paulo, Brazil, 2009. [29] X. Tong, X. Meng, K. Ding, Estimating geometric parameters of highways and railways using least-squares adjustment, Surv. Rev. 42 (318) (2010) 359–374. [30] S. Wolfram, Mathematica: A System for Doing Mathematics by Computer, Second ed. Addison-Wesley, Reading, MA, 1991. [31] L.L. Schumaker, Spline Functions: Basic Theory, Wiley-Interscience, New York, 1981. [32] P.M. Prenter, Splines and Variational Methods, A Wiley-Interscience Publication, New York, 1989. [33] C. de Boor, Spline Toolbox User's Guide, MathWorks, 2004. [34] A. Kouibia, M. Pasadas, Approximation by discrete variational splines, J. Comput. Appl. Math. 116 (1) (2000) 145–156. [35] A. Kouibia, M. Pasadas, A note on the discrete approximation of discontinuous curves and surfaces, J. Comput. Appl. Math. 208 (2) (2007) 373–379. [36] Ministerio de Fomento, Norma 3.1 IC. Trazado de la Instrucción de Carreteras, Ministerio de Fomento, 2000. [37] G. Koorey, Road data aggregation and sectioning considerations for crash analysis, Transp. Res. Rec. 2103 (2009) 61–68.