Integration of modern remote sensing technologies for faster utility mapping and data extraction

Construction and Building Materials xxx (2017) xxx–xxx

Contents lists available at ScienceDirect

Construction and Building Materials journal homepage: www.elsevier.com/locate/conbuildmat

Integration of modern remote sensing technologies for faster utility mapping and data extraction Aleksandar Ristic´ ⇑, Zˇeljko Bugarinovic´, Milan Vrtunski, Miro Govedarica, Dušan Petrovacˇki University of Novi Sad, Faculty of Technical Sciences, Department of Computing and Control Engineering, Trg Dositeja Obradovic´a 6, Novi Sad 21000, Serbia

h i g h l i g h t s
Integration of two NDT technologies, application on district heating network.
Extraction of coordinates of hyperbola apex and prongs.
Heating pipelines classification and geometric parameters estimation.
Tested on radargams with various pipes’ diameters.

a r t i c l e

i n f o

Article history: Received 3 November 2016 Received in revised form 19 May 2017 Accepted 3 July 2017 Available online xxxx Keywords: District heating network GPR Aerial thermography Neural networks Edge detection Automated data extraction

a b s t r a c t The aim of the research presented in this paper is to analyze the benefits of integrating a mobile system capable of very fast, reliable and relatively inexpensive detection, identification and status examination of district heating network. Thermal imaging using unmanned aerial vehicle is used for pipeline route detection, inspection of validity of cadastral data and for locating possible leakages. Ground Penetrating Radar – GPR technology is used for control sampling of radargrams on specific locations of routes in order to achieve following: identification of the geometric characteristics of district heating pipelines and structure, prevention and registration of damage, as well as automated data extraction. The main part of the paper is dedicated to the algorithm for automated data extraction, based on artificial neural networks and pattern recognition. Radargrams of district heating pipeline were used as input data for the extraction algorithm, while the results are geometric characteristics such as pipe depth, distance between pipes and diameter. Ó 2017 Elsevier Ltd. All rights reserved.

1. Introduction – Target characteristics District heating systems which transport hot water to consumers are often used in urban environments, because they are economically justified for large number of users [1]. In Serbia for instance, such systems are present in 57 cities and towns for more than 800,000 apartment units, which is 25% of the entire residential capacities, while in Novi Sad it is applied in more than 60% of apartment units [2]. The system consists of a heating power plant with additional objects, as well as of pipeline network. The main characteristic of district heating network is that the pipes are installed in pairs, mostly underground. Primary line is used to transport hot water (T > 110 °C) to the consumers, while the secondary line is used to transport water back to power plant, to be ⇑ Corresponding author. E-mail addresses: [email protected] (A. Ristic´), [email protected] (Zˇ. Bugarinovic´), [email protected] (M. Vrtunski), [email protected] (M. Govedarica), petrovacki@uns. ac.rs (D. Petrovacˇki).

heated again. Along with this pair of pipes, systems with three pipes are sometimes used, where the third pipe is used for distribution of hot sanitary water. Considering the fact that the transported medium is hot water, it is rather significant to take energy efficiency of the network into consideration, therefore a good quality thermal insulation of network is necessary. Older segments of network are constructed with seamless steel pipes installed in channels made of reinforced concrete with relatively small depth, covered with reinforced concrete panels, without the leakage alarm system (Fig. 1a). Channels are built on site, as reinforced concrete structures or out of prefabricated units with the dimensions ranging from 75  50 cm up to 250  150 cm. Pipes for heating lines are insulated later (independent insulation), most often with insulation materials based on glass mineral wool. Significant temperature changes of transported fluid influence the faster deterioration of pipes as well as significant dilatation changes of pipeline which are compensated by using construction compensator (‘‘U”, ‘‘L”, ‘‘Z” shape of compensator) or by some other design [2].

http://dx.doi.org/10.1016/j.conbuildmat.2017.07.030 0950-0618/Ó 2017 Elsevier Ltd. All rights reserved.

Please cite this article in press as: A. Ristic´ et al., Integration of modern remote sensing technologies for faster utility mapping and data extraction, Constr. Build. Mater. (2017), http://dx.doi.org/10.1016/j.conbuildmat.2017.07.030

2

A. Ristic´ et al. / Construction and Building Materials xxx (2017) xxx–xxx

Fig.1. Characteristics of district heating pipeline in concrete channel (a) and in trench (b).

The development of modern district heating systems showed that the usage of pre-insulated seamless steel pipes with leakage alarm system (Fig. 1b), put in an excavated trench and then covered with sand and dirt, represents a more reliable, more efficient solution in terms of energy and costs, as well as longer lifetime [1,3]. In many cities in Serbia, as well as in Novi Sad, there are reconstructions of the network underway. Old pipes are replaced with pre-insulated pipes, directly put into the trench and covered with sand. The diameter of the pipes which were used to convey hot water is determined by the importance of the distribution lines (primary, secondary and end-user connections). The geometric characteristics and dimensions of district heating pipelines in concrete channels and trenches are highly standardized [4], which allows facilitating the characteristic identification and the detection of routes. 2. NDT technologies for detection and analysis of district heating network Current technologies for non-invasive detection and analysis (NDT) of the state of district heating pipelines can be classified into technologies for detection and analysis of geometric characteristics of pipelines (Ground Penetrating Radar – GPR, thermal camera, electromagnetic locator – EML) and into technologies for pipeline state analysis (thermal camera, ultrasound methods, leakage alarm system, monitoring of flow and pressure within the pipeline) [5–7]. 2.1. NDT technologies for detection of district heating network Preparation of valid data for introduction of digital cadastre of utilities and its regular updating is one of major tasks for NDT technologies. Although there are many such technologies, dynamic development of urban areas requires technologies with fast acquisition, as well as very reliable identification of underground utilities. For instance, one of the recommendations for underground utilities detection in urban areas given by The Survey Association – TSA [8] is combined detection using EML and GPR technology. Considering specific geometry of heating pipeline, EML technology can be used for precise locating of pipeline routes and, in significant amount, determination of depth to the center of the pipe, which is the main characteristic of this technology [9]. Other geometric characteristics of heating pipelines presented in, as well as

the deterioration level of a pipe, cannot be determined with EML. Therefore, GPR technology, with the variety of its applications [10], represents the solution for identification of all geometric characteristics of a pipeline, especially when the documentation is incomplete, unreliable or even does not exist. The GPR is one of the most significant and advanced geophysical techniques [11] that detects changes in electromagnetic properties and provides high resolution data [12]. GPR can provide information on electromagnetic properties of the material as well as on underground objects. Electromagnetic properties can be found through fullwave inversion while high-resolution subsurface imaging can be used to locate objects [13]. In recent period, GPR technology has become more accessible and more present in engineering applications [10]. The result of data acquisition with GPR technology is radargram which contains the description, in visual, spatial and quantitative sense, of changes of transmitted signal within the time record of its reflection, after the propagation through heterogeneous medium containing various manmade objects. Those changes are caused by the existence of both metal and non-metal manmade objects (for instance, underground utilities), structures and voids within soil layers filled with air. Properties of surrounding soil, primarily dielectric permittivity, soil moisture content and soil homogeneity have impact on reflected signal as well. Determining correct dielectric constant of soil enables estimation of averaged wave propagation velocity and therefore calculation of relative depth of underground object. Using appropriate techniques of data processing it is possible to analyze the radargram in both time and frequency domain, and also perform identification of the cause of the change along with localization as well as estimation of their geometric and quality parameters [10]. If GPR acquisition is applied in zones with individual underground objects of cylindrical shape, a characteristic shape of hyperbolic reflection is formed in radargram, which provides an opportunity for development and application of automated methods for their detection in radargram [14]. Since in real field conditions congestion of underground utilities is often very high appearance of multiple hyperbolic reflections in radargram is expected. Furthermore, other interference such as clutter caused by spatial disposition of the underground utilities, heterogeneity of the surrounding medium and GPR device noise are also present. The complexity of hyperbola detection is also affected by asymmetric hyperbolic reflections, intersected, distorted, as well as incomplete hyperbolic reflections. Commercial software applications provide graphical characterization of the hyperbolic reflections along with wave velocity analysis [13].

Please cite this article in press as: A. Ristic´ et al., Integration of modern remote sensing technologies for faster utility mapping and data extraction, Constr. Build. Mater. (2017), http://dx.doi.org/10.1016/j.conbuildmat.2017.07.030

A. Ristic´ et al. / Construction and Building Materials xxx (2017) xxx–xxx

Existing procedures for radargram analysis can be performed by either analyzing the full, dense radargram image or by analyzing a thresholded sparse version of it [15]. It is possible to apply unsupervised procedures (Hough transform, for instance) and supervised procedures (e.g., Artificial Neural Networks – ANN) if dense radargram is analyzed [15]. Existing strategies for hyperbolic reflections detection involve application of algorithms that implement Hough transform [16–19], Wavelet transform [20], Radon transformation [21], standard algorithms for pattern recognition, such as Support Vector Machines (SVM) [22,23] or ANN [24]. The Hough transform is a computationally intensive method and has a cubic time complexity [15]. ANNs can be trained using signal processing statistical data descriptors [25], Welch power spectral density estimate of signal segments [26] or generated data sets [27,28]. Efficiency of ANN depends on characteristics and amount of training data [29]. In [13] the authors propose an automated method for identifying multiple hyperbolic reflections in inhomogeneous surrounding soil. The algorithm detects peaks of hyperbolae by setting analytically defined function of a hyperbola on a profile of edge points detected using the Canny filter. The existence of hyperbola is defined by a series of carefully selected criteria in order to match the zone of interest. A large number of methods is based on common characteristics of image segments, where the trained neural network is used [19,25,26,29–32]. The main objective of this approach is to narrow the search area. To reduce the need for a large database, the authors in [22] use simulated data for the training set. Also, in [24] the authors illustrate how different levels of noise can affect the performance of detection. The analysis of the papers related to these technologies yields the conclusion that the search through the dense radargram is very demanding in terms of time and computation resources and sensitive to noise and hyperbolic segments interference as well. Also, majority of the procedures analyze simplified radargrams. Simplification can be done in two ways [29]: simplification by using edge detection [33] or binarization [34] or segregation of small twodimensional sections from dense radargram (called segment of interest - SOI). Data from the hyperbolic reflections are then extracted from simplified radargrams. The objective of this paper is development of automated method for detection and point extraction from reflections of heating pipelines. Radargrams with reflections of heating pipes contain either two interfered hyperbolic reflections or reflection from reinforced concrete channel which is not hyperbolic. Therefore, methods mentioned above which can be used to detect multiple hyperbolic reflections are not completely suitable solution for case of heating pipelines. Proposed solution incorporates usage of trained ANN for automated detection of zones (SOI) containing specific shapes of reflections. After the search area is narrowed, high resolution image analysis is used to extract points from both types of reflections. The last step is estimation of geometric characteristics of district heating pipelines based on extracted points. In case of two pipes in trench estimation can be done by using fitting algorithm based on physical model of theoretical hyperbola [14], while in case of concrete channel estimation is done based on standardized dimensions of the channel and pipes. Proposed algorithm represents one possible complete solution for analysis of geometry of district heating pipelines. Algorithm was tested on a number of radargrams acquired in the field survey. To illustrate the performances of the proposed algorithm, we present the characteristics of the algorithm through seventeen representative radargrams. In these examples, we present different acquisition scenarios for both cases by varying pipe diameter, antenna frequency and acquisition parameters.

3

2.2. NDT for the state of district heating network analysis The term state of district heating system (the degree of pipeline deterioration) means the analysis of damages on insulation and the degree of corrosion of the pipes by determining the zones with damaged insulation and/or zones with leakages [35]. Since the degree of corrosion is difficult to estimate using NDT (especially during the time pipeline is in use), state of district heating system is determined by registering increased heat radiation at places with leakages or damaged insulation [36]. Large-scale methods for monitoring a district heating network are based on usage of aerial (airplane or unmanned vehicle) [37] and mobile terrestrial platforms with sensors (mounted on vehicle or hand-held devices). Aerial thermography is very convenient technology for the analysis of the state of district heating system. It involves the acquisition of images using unmanned sensor platform with thermal camera being the sensor. In comparison, satellite and airplane platforms have a drawback in terms of sensor resolution and costs. One of the advantages of unmanned platform is fast acquisition on quite large area of interest as well as lower costs, when compared to other technologies. General drawback of thermography is its short period of applicability during the year, which is when the heating system is operational. Aerial thermography has also the limitation in terms of registering some other objects with higher temperature than the surrounding (chimneys, buildings with heating losses, cars). Another constraint is that thermal image represents only current state of district heating pipeline, not the entire spatio-temporal change of the network [3]. Problem can be solved by developing a database of the state of district heating network and updating it in regular intervals. Updating intervals may vary from two years to taking images in the phase of testing and at the end of heating season (2 images per season). 3. Proposed solution description Based on the given analysis of available technologies, the combined application of two technologies selected in a way that they complement each other and overcome the drawbacks of those technologies taken separately arises as the most appropriate solution. The solution proposed in this paper involves the application of aerial thermography with unmanned vehicle for fast and efficient locating of heating pipelines as well as for the analysis of the state of district heating network, while examination of specific segments of the route and analysis of geometric characteristics of pipeline would be done using GPR technology. Locating heating pipelines by aerial thermography is based on the creation of georeferenced image where the color of each pixel corresponds to the temperature. Defining the position model of district heating pipeline network is done by selecting those pixels from the route which represent higher temperature compared to the pixels outside the route, with all the constraints mentioned above [5]. GPR technology is used for: a. Classification of district heating pipeline and estimation of geometric characteristics. When the database of district heating network is incomplete, unreliable or doesn’t exist a problem appears when pipeline has to be replaced or repaired. b. If heating pipeline is not operational then GPR technology represents appropriate way to identify, classify and estimation of geometry.

Please cite this article in press as: A. Ristic´ et al., Integration of modern remote sensing technologies for faster utility mapping and data extraction, Constr. Build. Mater. (2017), http://dx.doi.org/10.1016/j.conbuildmat.2017.07.030

4

A. Ristic´ et al. / Construction and Building Materials xxx (2017) xxx–xxx

c. Proposed procedure automates the process of identification of heating pipeline type and geometry, even when during the acquisition it is not possible to form the radargram with hyperbolic reflections (case of the concrete channel).

2. By localizing hyperbolic reflections, extracting coordinates that characterize the apices and the points on prongs of a hyperbolic reflection, as well as estimating parameter [14]. 3.1. Instrumentation

In Fig. 2 a standardized disposition of heating pipelines is given [4] in case of concrete channel (Fig.2a) and in case of trench (Fig. 2b), as well as measures that are later used to determine geometric characteristics. The identification of the type is done by examining the shape of reflection of heating pipes in radargram (concrete channel - Fig. 2c, trench - Fig. 2d). Estimated geometric characteristics of district heating network in this paper are diameter and a depth relative to the apex of the pipe, as well as dimensions of concrete channel. Since in a case of pipeline in concrete channel, due to concrete reinforcement, a high level of interference is present, determining pipelines dimensions is solved by estimation of concrete channels dimensions, and then, considering values in Table 1, pipes can be classified into three categories according to diameter (small, DN 65–80; medium, DN 100–150; large, DN 200–250) and depth is determined. It should be noted that, if there is a cover on a channel, the diameter is determined by estimating of concrete channels dimensions, and if cover is removed, the estimation based on distance between axes of the pipes is performed. If pre-insulated pipes in trench are to be examined, the problem can be solved in two ways: 1. By localizing the apexes of both hyperbolic reflections, estimating the axial distance which is standardized and determination of diameter according to values from Table 2.

The instrumentation used in those experiments included infrared thermographic camera mounted on unmanned vehicle and ground penetrating radar with two antennas of different central frequencies. Unmanned Aerial Vehicle Aibot X6 is used as a platform for thermal camera. It was controlled remotely from the ground. The vehicle has the flight time of 30 min. Mountable sensors are, besides thermal camera, digital camera, multispectral camera etc. The Optris PI Lightweight Kit consists of lightweight PC and infrared camera PI450 in lightweight housing. Kit is developed specifically for use in aerial thermography. All measurements are created with a time-domain radar, the model SIR3000 from Geophysical Survey Systems, Inc. (GSSI, Salem, MA, USA). Two shielded bowtie antennas were used, with 200 MHz (Model 5106) and 900 MHz (Model 3101) center frequency. Scanning parameters are given in Section 6. Control unit and antennas were mounted on a vehicle equipped with incremental encoder so all GPR scanning is done in distance mode. 3.2. Contribution One of the main goals of the research presented in this paper is the proposal of integrated technology for efficient acquisition and monitoring of the state of district heating network, based on

Fig. 2. Geometry of heating pipes in concrete channel (a) and in trench (b). Radargram done with 200 MHz antenna containing reflections of heating pipeline in concrete channel (c); Radargram done with 900 MHz antenna containing reflections of heating pipeline in trench (d).

Please cite this article in press as: A. Ristic´ et al., Integration of modern remote sensing technologies for faster utility mapping and data extraction, Constr. Build. Mater. (2017), http://dx.doi.org/10.1016/j.conbuildmat.2017.07.030

A. Ristic´ et al. / Construction and Building Materials xxx (2017) xxx–xxx

5

Table 1 Dimensions of concrete channel [4]. DN

b [cm]

Bi [cm]

Bc [cm]

Dp [cm]

Hi [cm]

Hc [cm]

65 80 100 125 150 200 250

20 20 23 23 23 55 55

65 65 80 80 80 125 125

85 85 100 100 100 145 145

76 89 108 133 159 216 267

50 50 50 50 50 65 65

70 70 70 70 70 85 85

Table 2 Dimensions of trench [4]. DN

D x s/Dp [cm  cm/cm]

B min [cm]

b min [cm]

Bu min [cm]

H min [cm]

Hu min [cm]

65 80 100 125 150 200 250

7.61  0.29/14 8.89  0.32/16 10.8  0.36/20 13.3  0.40/22.5 15.9  0.45/25 21.91  0.63/31.5 27.3  0.63/40

15

29 31 35 37.5 40 56.5 65

63 67 75 80 85 108 125

50

74 76 80 82.5 85 91.5 100

25

analysis of current state in research area and performed experiments. Integrated technology is a combination of two remote sensing technologies: thermal camera mounted on unmanned aerial vehicle and GPR technology. The second task, to which the most of the attention in this paper was dedicated, is the development of procedure for automated identification and classification of district heating network systems, as well as for estimation of geometric characteristics of pipelines in terms of determining pipes’ diameter and depth relative to the pipes’ apex. In proposed algorithm, training of ANN was done by using the existing Matlab function CascadeObjectDetector with our original training set. Extraction of points of the apex and on the prongs of hyperbolic reflection/concrete cover is original contribution of this paper. The procedure represents one possible solution for automated analysis of district heating pipelines.

4. Artificial neural networks For a successful data extraction trained neural network was used to analyze radargram in the form of raster image, resulting with extracted segments of interest (SOI). This significantly reduces the amount of data for further analysis. Implementation of automated SOI selection in this paper was done by using the existing function CascadeObjectDetector (COD) which is the part of the Computer Vision System Toolbox – Object Detection and Recognition in Matlab R2014a. COD can detect objects whose aspect ratio is not significantly different. It is based on machine learning and it demands a training network on its input. In case of supervised classification, as in this paper, a set of input data is created based on selected positive training samples containing object of interest, and on selected negative training samples, not containing object of interest. Quality control of a trained classifier is done based on input test data. If output values do not solve entirely a problem previously defined, it is necessary to use larger data sets, in order to provide better detection ability so that the classifier can generalize new data sets better. Due to different geometry of the object of interest, in cases of district heating pipes in concrete channel and in trench, two training networks were created. The first one is made using samples representing pipes in concrete channel and the second one using

samples representing pipes in trench. The algorithm is searching through the input raster image by moving the search window across its smaller regions. After that, the detector uses cascade classifier to decide whether the current window contains the object of interest or not. The size of the window is changing so that the detector can recognize the same object in different scales, but it is important that the aspect ratio of object of interest remains the same. COD consists of several phases and in each phase the region of the search window is marked as positive or negative. Phases are designed in such manner that negative regions are discarded as soon as possible, presuming that the most of the windows do not contain the object of interest. Search windows marked as positive in one search phase are forwarding the possible positive solution to the next phase; otherwise the detector moves the window to next location. One of following models can be used for training within COD:
Haar-like features.
LBP – Local Binary Patterns.
HOG – Histogram of Oriented Gradients. Haar model utilizes the principle of the analysis of a number of pre-defined frames when detecting objects and divides them in two classes: object and non-object. It is most frequently used model for face detection [31]. When LBP model is used, central pixel compares its value with values of 8 neighboring pixels, in circular order, clock- or counter clockwise. If the value of neighboring pixel is higher, 1 is entered to that position, otherwise 0. That way a binary eight-digit number is created and then it is converted into decimal one due to easier comparison [38]. HOG is based on the fact that occurrence of local object and its shape within the image can be described using the ordering of intensity gradient, or by directions of edges. The original image is divided into smaller connected regions (cells) and for each cell either a histogram of gradient directions is calculated or the orientation of edges for pixels that are within the cell. To improve the performance, the contrasts of local histograms can be normalized by calculating the intensity over larger region of image (block) and then by using that value to normalize all cells within the block [39]. After the analysis of all three models, HOG turned out to be most successful solution and therefore is used to obtain the results represented in this paper. Quality criterion for each model was the rate of false detected objects. Object is considered to be false

Please cite this article in press as: A. Ristic´ et al., Integration of modern remote sensing technologies for faster utility mapping and data extraction, Constr. Build. Mater. (2017), http://dx.doi.org/10.1016/j.conbuildmat.2017.07.030

A. Ristic´ et al. / Construction and Building Materials xxx (2017) xxx–xxx

6

positive when negative region is misclassified as positive, and false negative when positive region is misclassified as negative. Overall s false negative rate of cascade classifier is f , where f is false alarm rate per phase in range (0,1) and s number of phases. Similar, overall true positive rate is ts, where t is true positive rate per phase in range of (0,1]. Therefore, it can be noticed that adding phases means lower overall false positive rate, but decreases overall true positive rate [40]. Computation time for Haar model was slightly above 2 min, while for LBP and HOG model it was around 30 s or less. Number of detected frames for approximate locating of objects was 1 for HOG model, while for LBP and Haar model it was 13 and 14, respectively (Fig. 3). For creation of training network, the parameters were adjusted to yield optimal results. They are mutually dependent, so increasing the number of phases decreases the false detection rate but increases the probability for some reflection not to be found. Table 3 contains parameters used in this training. Emphasis in this paper is not on extraction of SOI, but on extraction of geometry data of observed object in SOI.

5. Algorithm description The algorithm developed for detection of objects of interest (OOI) and estimation of parameters specific for heating pipelines consists of several steps (Fig. 4). Different shapes of reflection in radargram for pipes in concrete channel and in trench indicate different methods for data extraction. In the first step of the algorithm, approximate locating of OOI is performed, using training network and decreasing the size of input raster image to smaller SOI which may contain searched object. After this step, it is necessary to extract characteristic data while minimally relying on resulting boundaries of the frame. Considering high level of interference, which is the consequence of geometry and reinforcement of concrete channel, appearance of hyperbolic reflections cannot always be expected. Therefore, the focus of the algorithm is directed towards extraction of points which represent concrete channel cover, aiming to enable the assessment of the width of the cover, that is the channel itself. Comparing assessed width with the values in Table 1 it is possible to determine the group of pipes diameters for considered channel. For instance, if the cover width is 100 cm, according to Table 1 possible diameters are DN100, DN125 and DN150. Determination of exact diameter (like in the procedure for the pipes in trench) is

possible only if there is at least one usable hyperbolic reflection of the pipe. In case of preinsulated pipes in the trench, results of the application of this algorithm are extracted point from hyperbolic reflections, with defined coordinates of both apices (axial distance of two pipes). These results represent input data for algorithm described in [14] where which simultaneously estimate cylindrical object radius and averaged wave propagation velocity. Since the axial distance is determined using proposed algorithm, applying standardized data from Table 2 estimated radius of the pipe can be checked. Estimated value of averaged wave propagation velocity is used to determine the depth to the pipes apex. 5.1. Heating pipeline in concrete channel – Analysis in time domain In case of pipeline in concrete channel, bounding frame of selected SOI is used to start of the search across the width of the channel. Large part of EM waves is reflected of the concrete channel so the pipeline has to be determined indirectly via the geometry of the concrete channel while in case of the trench the detection is done directly using characteristic hyperbolic reflections. The width of concrete channel is determined over the cover plate because that is where the reflection is most noticeable. The first step of the search is finding a pixel with the maximum intensity within the middle column of SOI (Fig. 5). Considering that in lower part of SOI quite often pixels with high intensity can be found, algorithm finishes the search when it finds the first peak of the column which is more noticeable and which is accepted as the starting pixel of the search. For the first iteration, search of local maximums, the algorithm extracts the part of input raster to the left side so that the starting pixel is in second row of the second column in submatrix 3  2. After that the search is moving to both left and right side by finding local maximums within 3  2 submatrix, without current pixel with maximum intensity. Submatrix is moving iteratively so that the current maximum is always in second row of the second column of extracted submatrix. Search on the right side has identical approach for finding pixels, but current maximum is always in second row of first column of extracted submatrix (Fig. 5). T is the current pixel with maximum intensity. In order to determine the border of area where pixels with local maximum are being searched for, two constraints are adopted. The first one is related to recognizing the edge points of concrete cover plate based on movement of the submatrix. Summing is performed when the submatrix, while going from one iteration to the next

Fig. 3. Resulst of Haar model (a), LBP model (b) and HOG model (c).

Table 3 Optimal training parameters. District heating network

False alarm rate

True positive alarm rate

Number of stages

Model

Number of negative samples

Concrete channel Trench

0.44 0.40

0.995 0.995

10 9

HOG HOG

4 1

Please cite this article in press as: A. Ristic´ et al., Integration of modern remote sensing technologies for faster utility mapping and data extraction, Constr. Build. Mater. (2017), http://dx.doi.org/10.1016/j.conbuildmat.2017.07.030

A. Ristic´ et al. / Construction and Building Materials xxx (2017) xxx–xxx

7

Fig. 4. Basic steps of the search algorithm and extraction of data characteristic for heating pipelines in concrete channel and in trench.

Fig. 5. Determination of concrete channel width.

one, is moving downwards or remains in the same row. It is being performed until the submatrix moves upwards (and cancels itself to zero) or reaches boundary value. After the analysis of large number of radargrams containing the reflection of the concrete channel with known dimensions, boundary value of mentioned sum, is determined to be 4. Therefore, when algorithm reaches provided value the further search is aborted.

For the second constraint, the value of pixel intensity is used. The value of current local maximum is compared to the value of intensity of starting pixel. If their ratio is lower than 87%, the further search is aborted. Value of 87% is determined experimentally, applying proposed algorithm on a number of radargrams containing reflections of concrete covers of known dimensions. Increasing this value would yield detected dimensions smaller than real ones,

Please cite this article in press as: A. Ristic´ et al., Integration of modern remote sensing technologies for faster utility mapping and data extraction, Constr. Build. Mater. (2017), http://dx.doi.org/10.1016/j.conbuildmat.2017.07.030

A. Ristic´ et al. / Construction and Building Materials xxx (2017) xxx–xxx

8

while decreasing would yield bigger dimensions. Whichever constraint is violated further search is aborted. Extracted pixels are then used to determine the dimensions of the channel (width and height) and the depth to the cover. Indexes of columns of the pixels on the far right (Cr) and far left side (Cl) are subtracted and then divided with the value of the scanning parameter scan/m resulting in width (W) (Eq. (1)). According to determined width and values from Table 1 height of the channel can be assessed.

W¼

Cr  Cl ½m scans=m

ð1Þ

To determine the depth of the cover (D), total depth (Dtot) is calculated as in Eq. (2).

c Range ½m Dtot ¼ pffiffiffiffi  er 2

ð2Þ

Also, average value of row indexes of extracted pixels is calculated (Rav). Total depth is divided by numbers of samples per scan to calculate the value of sample in meters. Those two values are then multiplied resulting with the value depth (Eq. (3)).

D ¼ Rav 

Dtot ½m Samples=scan

Fig. 6. Magnitude image DN250.

ð3Þ

5.2. Heating pipeline in concrete channel – Frequency domain analysis For a more complete analysis of radargrams that contain reflection from the concrete channel, frequency domain was also taken into consideration. Conversion from time to frequency domain is done using Fourier transformation. Most frequently, data are separated into time segments of a certain duration, and the transformation is subsequently done for each segment. Optimum size of time window for radargram analysis can be calculated directly from known scanning parameters (Eq. (4)). Number of points in time for which the spectrum is calculated depends on the size of the window and the degree of overlap between neighboring windows. The overlap of windows causes the correlation of neighboring results. As suggested in [40], recommended value for overlap is 50%. Spectrums for all time windows form a time-frequency representation of a signal for each scan called spectrogram [40].

win ¼ 2  T ant ¼

2b Af  range

Fig. 7. Frequency imageDN250.

ð4Þ

win – Optimum size of time window [samples]. b – number of samples per scan. Af – antenna central frequency [MHz]. range – two-way travel time [ns]. Magnitude image can be created based on spectrogram. It is important for forming binary mask and extracting the zone of interest from the radargram. Each point on magnitude image represents the value of maximum magnitude of scan spectrogram at a certain location in a certain moment [40]. The magnitude image clearly shows that the reflections of the highest intensity originate from the parts of a radargram that contain the reflections of heating pipelines. Combining magnitude and frequency image (Figs. 6 and 7) a mask is obtained for extracting segments of interest from the radargram. The mask is of a binary type, with the value 1 assigned to the segments of radargram that meet two conditions: 1. Values of magnitude within the range of 80–100% compared to maximum value from the magnitude image. 2. Frequencies in the range of 120–220 MHz (frequency image).

The resulting binary mask (Fig. 8) applied on radargram by element-wise multiplication extracts segments of radargram which may contain object of interest (Figs. 9a and 10a). After the extraction, Canny algorithm for edge detection is applied on obtained segments. It finds the edges by local maximum of image (radargram) gradient, which is calculated using Gauss filter. A filter mask is used and the larger the mask, the lower is detector’s sensitivity to noise. On the other hand, increasing the mask yields to bigger error of edge localization. Radargram gradient is calculated after filtering. The method uses two threshold values in order to detect strong and weak edges, and includes weak edges into output result only if they are connected with strong edges. The result of the detection in comparison to the real radargram is shown in Figs. 9b and 10b. The Matlab function houghpeaks (Image Processing Toolbox) is used to define rows and columns derived after Hough transformation. This function sorts newly obtained pixels and assigns them values for x and y coordinates. After that, forming lines from segments (pixels) is done using function houghlines (Image Processing Toolbox). Figs. 9c and 10c show the processing results of the radargram for the pipe DN250 and DN80 respectively. The final form of binary radargram obtained through the described procedure is

Please cite this article in press as: A. Ristic´ et al., Integration of modern remote sensing technologies for faster utility mapping and data extraction, Constr. Build. Mater. (2017), http://dx.doi.org/10.1016/j.conbuildmat.2017.07.030

A. Ristic´ et al. / Construction and Building Materials xxx (2017) xxx–xxx

9

yielding the value of a sample in centimeters. This value is then multiplied with the difference in samples between horizontal lines corresponding to top and bottom of the concrete channel to obtain its height. 5.3. District heating pipelines in trench

Fig. 8. Binary mask.

used as a background in figures. The width and the height of concrete channel can be determined based on extracted lines. Width is calculated by dividing the length of the line in pixels with the value of the scanning parameter scans/meter. To determine the height, velocity has to be known. Range is divided by a number of samples and resulting value [ns] is multiplied with velocity [ns/cm]

In order to calculate the elements that characterize the pipeline in a trench it is necessary to determine the apexes of hyperbolic reflections for two pipes. The procedure is performed in two steps. First, coordinates of apices of hyperbolic reflections are extracted and then points on prongs (hyperbola branches). Apices are used as starting points for extraction of points on prongs. In the first step of the search, extracted segment is divided into four equal parts and its boundaries are used subsequently (Fig. 11). The initial pixel to start the search is found in column between the first and second quarter of SOI, and second one between third and fourth quarter of SOI. The coordinates of pixels with maximum values of intensity in previous columns are considered to be the coordinates of the initial pixels. After that, 3  2 submatrices are formed in the same manner as for pipes in concrete channel. Since the pipes are close to each other, the prongs of hyperbolic reflections cross. Having that in mind, the search has to be stopped before the crossing point, in order to prevent the search to go to the next hyperbolic reflection. Therefore, the stopping criteria is defined based on submatrix movement, with the limit value for the sum to be six. This way, pixels that are extracted represent plausible apex of hyperbolic reflection. After the search stops it is necessary to adopt unique

Fig. 9. Radargram containing pipeline DN250 (a), Canny edge detection (b), Hough function (c).

Please cite this article in press as: A. Ristic´ et al., Integration of modern remote sensing technologies for faster utility mapping and data extraction, Constr. Build. Mater. (2017), http://dx.doi.org/10.1016/j.conbuildmat.2017.07.030

10

A. Ristic´ et al. / Construction and Building Materials xxx (2017) xxx–xxx

Fig. 10. Radargram containing pipeline DN80 (a), Canny edge detection (b), Hough function (c).

Fig.11. Determination of unique apex of hyperbolic reflection.

coordinates of the apex. For the unique row, the one with the lowest index of previously extracted pixels is accepted. To determine the unique column of the apex, extracted pixels are divided into two groups (left and right prong), based on determined row. Pixels that are three or four rows away from the accepted row are found. If there is more pixels in the same row within one group their

mean value of index is calculated, determining coordinates of pixel 1, 2, 3 and 4 (Fig. 11). Later on, the algorithm finds the values of approximate columns Sr1 and Sr2, and as unique column of the apex their mean value is accepted (Eq. (5)). The same procedure is applied for both hyperbolic reflections, with the opposite directions of search.

Please cite this article in press as: A. Ristic´ et al., Integration of modern remote sensing technologies for faster utility mapping and data extraction, Constr. Build. Mater. (2017), http://dx.doi.org/10.1016/j.conbuildmat.2017.07.030

A. Ristic´ et al. / Construction and Building Materials xxx (2017) xxx–xxx

11

Fig. 12. Extraction of points from hyperbolic reflection prongs.

Fig. 13. Test area on orthophoto image (a) and cadastral map (b).

Fig. 14. Location 1.

ApexColumn ¼

Sr1 þ Sr2 2

ð5Þ

The depth and distance between pipes’ axes are calculated after unique coordinates of apexes are found. Some parameters that were set when GPR scanning was done, are used in this part of algorithm, such as scans/m, sample/scan, range and dielectric con-

stant. Assessment of the depth of pipes in trench is similar to procedure applied to concrete channel in time domain analysis. The only difference is that here row index of apex is used (Rapex), so there is no need to calculate the average value (Eq. (6)).

D ¼ Rapex 

Dtot ½m Samples=scan

ð6Þ

Please cite this article in press as: A. Ristic´ et al., Integration of modern remote sensing technologies for faster utility mapping and data extraction, Constr. Build. Mater. (2017), http://dx.doi.org/10.1016/j.conbuildmat.2017.07.030

12

A. Ristic´ et al. / Construction and Building Materials xxx (2017) xxx–xxx

Fig. 15. Locations 3 and 4.

Fig. 18. Radargram ‘‘FILE_003”, diameter 2xDN125.

Fig. 16. Locations 1 and 5.

Fig. 19. Radargram ‘‘FILE_157”, diameter 2xDN250.

Fig. 17. Radargram ‘‘FILE_161”, diameter 2xDN80.

The last step of the algorithm extracts the coordinates of points on prongs of hyperbolic reflections. The search is similar as for determining the apex coordinates, with search window dimensions 2  2. The extraction starts from the apex, and prongs are analyzed separately. For the left prong 2  2 submatrix is formed so that apex in first iteration is in the first row and second column. Among other pixels local maximum is found which is used to move submatrix in next iteration so that current maximum is always in the first row of the second column.

Fig. 20. Radargram ‘‘FILE_171”, diameter 2xDN250.

Please cite this article in press as: A. Ristic´ et al., Integration of modern remote sensing technologies for faster utility mapping and data extraction, Constr. Build. Mater. (2017), http://dx.doi.org/10.1016/j.conbuildmat.2017.07.030

A. Ristic´ et al. / Construction and Building Materials xxx (2017) xxx–xxx

13

Fig. 23. Radargram ‘‘FILE_317”, diameter 2xDN100.

Fig. 21. Radargram ‘‘FILE_268”, diameter 2xDN250.

Fig. 24. Radargram ‘‘FILE_319”, diameter 2xDN65. Fig. 22. Radargram ‘‘FILE_295”, diameter 2xDN250. GPR measurements were done afterwards on five locations, using 200 MHz antenna:

The coordinates of the points on right prong are found identically, provided that the current maximum is always in the first row and first column of the submatrix (Fig. 12). This way provides successful extraction of points on prongs even beyond crossing point of hyperbolic reflection prongs (Fig. 12 – right side). 6. Experimental results 6.1. Thermal camera In order to examine the capabilities of integrate system, a part of district heating pipeline route (orange lines) was recorded using thermal camera mounted on Unmanned Aerial Vehicle (UAV). An area within university campus in Novi Sad was chosen as a test area and cadastral plans in scale 1:500 are examined before the test to determine the supposed pipeline routes (Fig. 13). Within 15 min of the planned flight routes of district heating pipelines in the area of 50000 m2 were recorded (images where corresponding temperature value is assigned to each pixel).

& & & & &

Location Location Location Location Location

1 2 3 4 5

– – – – –

DN250mm. DN80mm. DN250mm. DN65mm. DN250mm.

Those locations where examined more thoroughly on thermographic images. District heating pipeline route in the vicinity of L1 is given in Fig. 14. The dark stripe on location is from a car parked above the pipe. Intensive heat radiation can be noticed and the soil above heating pipeline is 3–4 °C warmer than the surrounding soil. One significant difference between real and cadastral data is present - no pipeline DN50. Heating pipeline route in the vicinity of L3 and L4 is given in Fig. 15. Different pipe diameters are present on those locations, so heat radiation is evident only on the pipe of the biggest diameter (DN250). In this image just as in previous one, the difference between real and cadastral data may be noticed - pipeline DN65 route does not match the one in the cadastral map. Also, intensive radiation (8–9 °C more than surrounding soil) indicates damage of insulation and possible leakage. Heating pipeline route with locations L1 and L5 is given in Fig. 16. The area in the focus of the camera appears to be slightly hotter, which may indicate some energy losses caused by insulation deterioration.

Please cite this article in press as: A. Ristic´ et al., Integration of modern remote sensing technologies for faster utility mapping and data extraction, Constr. Build. Mater. (2017), http://dx.doi.org/10.1016/j.conbuildmat.2017.07.030

A. Ristic´ et al. / Construction and Building Materials xxx (2017) xxx–xxx

14

Table 4 The result for the district heating system in a concrete channel - Time domain. Radargram Name

True width (T) [cm]

Antenna frequency [MHz]

Time domain The resulting width (M) [cm]

Max difference (T-M) [cm]

Medium depth cover [cm]

Medium processing time [s]

FILE_161 FILE_003 FILE_157 FILE_171

85 85 145 145

200 200 200 200

90 88 152 136

5 3 7 9

36 73 18 36

0.089 0.165 0.093 0.092

Table 5 The result for the district heating system in a concrete channel - Frequency domain. Radargram Name

FILE_161 FILE_157

True width (T) [cm]

Antenna frequency [MHz]

85 145

200 200

Frequency domain The resulting width (M) [cm]

Max difference (T-M) [cm]

Medium processing time [s]

88 144

3 1

2.86 2.42

Table 6 The result for the district heating system in a trench. Radargram name

True diameter pipe [mm]

Antenna frequency [MHz]

The resulting depth of the pipe - left [cm]

The resulting depth of the pipe - right [cm]

Distance between axes [cm]

Medium processing time [s]

FILE_268 FILE_295 FILE_317 FILE_319

DN250 DN250 DN100 DN65

200 900 900 900

87 89 73 74

86 87 73 73

67 65 40 35

0.134 0.109 0.075 0.084

Table 7 Results of testing with different scanning parameters. Radargram name

scans/ m

samples/ scan

bits/ sample

FILE_341 FILE_340 FILE_329 FILE_344 FILE_343 FILE_335 FILE_334 FILE_336 FILE_337

50 50 50 100 100 100 100 100 100

256 512 512 256 256 512 512 1024 2048

8 8 16 8 16 8 16 8 8

Depth

Depth error

Apex 1 [cm]

Apex 2 [cm]

Apex 1 [cm]

Apex 2 [cm]

88.1 88.6 86.9 89.0 89.0 85.9 86.4 88.2 88.0

90.0 90.5 88.3 90.5 90.9 85.9 87.6 88.9 87.9

+0.1 +0.6 1.1 +1.0 +1.0 2.1 1.6 +0.2 +0.0

+2.0 +2.5 +0.3 +2.5 +2.9 2.1 0.4 +0.9 0.1

6.2. Algorithm testing – Radargrams for concrete channel and trench Experimental results are shown through several representative radargrams for two types of heating pipelines analyzed in the paper, pipelines in a concrete channel and in trench. In case of pipes in a concrete channel acquisition is performed with the antenna of central frequency of 200 MHz, and for pipes in trench with antennas of 200 and 900 MHz. Analysis was done for smaller (DN65 - Fig. 24, DN80 - Fig. 17), medium (DN100 - Fig. 23, DN125 - Fig. 18) and larger diameters of heating pipelines (DN250 - Figs. 19–22). The algorithm has been successfully estimated the width of the concrete channel as in time and in frequency domain. The maximum difference in widths of the concrete channels obtained in the time domain was 9 cm. Average depth of the concrete-reinforced panels in most cases was up to 36 cm, except for the case of a channel that runs under the road (Radargram ‘FILE_003’), which was 73 cm. All analyses in the paper were carried out on a PC with Intel Core i3-4130 CPU with 3.40 GHz with MATLAB environment. Radargrams processing time in the time domain was below 0.2 s. In Table 4 it can be noticed that the smaller absolute difference between the actual and the estimated width of the channel is obtained for smaller diameters, while the maximum difference for larger channel dimension was 9 cm. Width of the channel in radargram ‘FILE_003’ was also measured in the field and it was 85 cm. Although the error at channels of larger dimensions is bigger, it does not

Axial distance [cm]

Axial distance error [cm]

Processing time [s]

74 72 66 73 71 66 65 69 69

+6 +4 2 +5 +3 2 3 +1 +1

0.051 0.072 0.057 0.054 0.086 0.069 0.085 0.099 0.143

prevent the conclusion that the width of the channel is 145 cm. It should be noted that the parameter scans/m during the acquisition was set to the value of 50, but the error of 1 scan enters a wrong assessment of the width of the concrete channel of 2 cm. Analysis of heating pipelines placed in concrete channel was done in frequency domain as well (Table 5). Estimated width of the channel was different from real value for only 3 cm. This leads to the conclusion that the value of error provides unambiguous classification of detected heating pipelines. In this case, 30 mm of error occurred at the pipelines of smaller diameter, in concrete channels 80 cm wide. For channels with 145 cm of width, error was 1 cm. Average processing time for frequency domain analysis was larger than for time domain, and was around 3 s. Heating pipelines DN250 in trench were scanned using antennas with frequencies 200 and 900 MHz. Since it is a sandy soil, pipes were successfully detected with both antennas. High intensity reflection can be seen in radargram acquired with antenna of 900 MHz, primarily due to higher level of radial resolution with respect to the antenna of 200 MHz. In both cases the algorithm successfully determined the apexes and the coordinates of points on the prongs of hyperbolic reflections. Also, the resulting value of axial distance is almost identical, differing by 2 cm. The resulting axial distance corresponds to values in Table 2 yielding the conclusion that those are larger diameter pipes. The depth of the left and right apexes is larger for 2 and 3 cm respectively, with radargrams aqcquired with antenna of 200 MHz (Fig. 21). Crossing zone of the adjacent hyperbolic reflections was more clearly determined to radargram acquired with 900 MHz antenna (Fig. 22).

Please cite this article in press as: A. Ristic´ et al., Integration of modern remote sensing technologies for faster utility mapping and data extraction, Constr. Build. Mater. (2017), http://dx.doi.org/10.1016/j.conbuildmat.2017.07.030

A. Ristic´ et al. / Construction and Building Materials xxx (2017) xxx–xxx District heating pipelines of smaller diameters laid in trench were scanned using antenna of 900 MHz, primarily due to higher level of radial resolution, and weak reflection in case the scanning the pipes DN250 using antenna of 200 MHz. In the case of pipes with diameter DN100 estimated depths of the left and right apexes of detected pipes are identical, while for pipes of diameter DN65 left apex is slightly below the right apex (Table 6). For heating pipelines of smaller diameters, differences between axial distances are smaller, so incorrect determination of diameter is more probable. For diameter DN65 determined axial distance was the same as for DN100 so it could be concluded that the pipes are of medium diameter. Axial distance determined for pipes of diameter DN100 was just as the limit value for DN150, so it lead to correct conclusion that the pipes are of medium diameter. For heating pipes of larger diameters differences in axial distances are bigger and it provides the estimation of real diameter of the pipes, as for pipes DN250 (Table 6). 6.3. Algorithm testing – Change of acquisition parameters Two heating pipes DN250 were scanned with a 900 MHz antenna. Pipe axes were at a distance of 68 cm and both pipes were buried in the trench at a depth of 88 cm. The pipe depth and the axial distance were measured when the trench was open. Two interfered hyperbolic reflections can be observed in the radargrams (Fig.2d). The following acquisition settings were varied, which affect the horizontal and vertical radargram resolution: [scans/m], [samples/scan] and [bits/sample] (Table 7). The time window was always 23 ns. Quantitative results of the experiment are presented in Table 7. For all input data, the algorithm successfully determined the coordinates of apices of both hyperbolic reflections as well as the points on their prongs. Different settings yielded minor differences in the obtained values of pipe depth and distance between hyperbola apices. From Table 7 it can be seen that the error for the axial distance was in the range 3 cm to +6 cm. Axial distance estimation error depends on horizontal scanning resolution [scans/m], but indirectly from vertical scanning resolution [samples/ scan] as well. Therefore, when data acquisition is done in order to provide input data for proposed solution scanning parameter values that can yield best results are 100scans/m and vertical resolution of 512samples/scan or higher. The error of depth estimation is within range from 2.1 to +2.9 cm. Although this error is acceptable, higher values of both vertical scanning resolution [samples/scan] and number of bits yield lower error of estimated depth.

7. Conclusion In this paper, an innovative approach to analysis of district heating networks is presented. It is a solution that integrates two NDT technologies, combining their advantages and overcoming their drawbacks. Analysis of the current NDT technologies was done, and after taking into consideration characteristics of district heating networks, two technologies were selected: thermal imaging using unmanned aerial vehicle and ground penetrating radar. The first one enables fast acquisition, route examination and leakage detection, while the other one enables estimation of geometric characteristics of a pipeline. Integrated solution is used to obtain data for experimental verification. The second part of the paper is related to automated procedure for identification and classification of heating pipeline as well estimation of geometry characteristics i.e. the diameter and depth. Procedure is applicable to both types of heating pipeline installations and therefore represents a universal solution in its domain. Experimental verification of the procedure is performed on radargrams acquired in real conditions and shown in paper with several representative ones containing reflections from the pipes of different diameters. Regarding GPR technology, analysis of district heating network represents a complex problem which has not been treated in significant degree so far. Therefore, the represented procedure is a unique solution for automated analysis of district heating pipelines. Acknowledgments This work has benefited from the network activities carried out within the EU (European Union) funded COST (European Coopera-

15

tion in Science and Society) Action TU1208 ‘‘Civil Engineering Applications of Ground Penetrating Radar”. References [1] C. Delmastro, E. Lavagno, L. Schranz, Underground urbanism: master plans and sectorial plans, Tunn. Undergr. Space Technol. 55 (2016) 103–111. [2] Serbian Chamber of Engineers, Thematic chapter 9.2: District heating and preparation of sanitary hot water, http://www.ingkomora.org.rs/eng/?id= indexeng, (accessed July 2016). [3] J. Rasmussen, J. Morsing, A.Dyrelund, How to ensure a long lifetime of the district heating networks, Hot & Cool magazine no.1/2016, http://www. senerges.rs/images/Zivotni%20vek.pdf (accessed July 2016). [4] Novi Sad DC PUC, District heating network of Novi Sad: Technical system, http://nstoplana.rs/tehnicki-sistem/, (accessed June 2016). [5] A. Berg, J. Ahlberg, M. Felsberg, Enhanced analysis of thermographic images for monitoring of district heat pipe networks, Pattern Recogn. Lett. (2016), http:// dx.doi.org/10.1016/j.patrec.2016.07.0 (in press). [6] A. Ristic, M. Govedarica, M. Vrtunski, D. Petrovacˇki, Integration of modern remote sensing technologies for faster utility mapping and data extraction, 2015, European Geosciences Union General Assembly 2015, GI3.1 Session ‘‘Civil Engineering Applications of Ground Penetrating Radar”, Vol. 17., EGU2015-2365-1, Vienna, Austria, April 12th - 17th, 2015. [7] S.B. Costello, D.N. Chapman, C.D.F. Rogers, N. Metje, Underground asset location and condition assessment technologies, Tunn. Undergr. Space Technol. 22 (2007) 524–542. [8] The Survey Association, The Essential Guide to Utility Surveyors, Issue 4, 2011, http://www.tsa-uk.org.uk (accessed July 2016). [9] British standard PAS128, Specification for underground detection, verification and location (2014) [10] A. Benedetto, L. Pajewski (Eds.), Civil Engineering Applications of Ground Penetrating Radar, Springer, 2015. [11] D.J. Daniels (Ed.), Ground Penetrating Radar, second edition., The Institution of Electrical Engineers, London, United Kingdom, 2004. [12] H.M. Jol (Ed.), Ground Penetrating Radar: Theory and Applications, Elsevier, 2009. [13] L. Mertens, P. Persico, L. Matera, S. Lambot, Automated detection of reflection hyperbolas in complex GPR images with No A priori knowleadge on the medium, IEEE Trans. Geosci. Remote Sens. 54 (1) (2015) 580–596. [14] A. Ristic, D. Petrovacki, M. Govedarica, A new method to simultaneously estimate the radius of a cylindrical object and the wave propagation velocity from GPR data, Comput. Geosci. 35 (8) (2009) 1620–1630. [15] R. Janning, A. Busche, T. Horváth, L. Schmidt-Thieme, Buried pipe localization usingan iterative geometric clustering on GPR data, Artif. Intell. Rev. 42 (3) (2014) 403–425. [16] C.G. Windsor, L. Capineri, P. Falorni, The Estimation of Buried Pipe Diameters by Generalized Hough Transform of Radar Data, Progress In Electromagnetics Research Symposium, Hangzhou, China, 2005, pp. 345–349 [17] G. Borgioli, L. Capineri, P. Falorni, S. Matucci, C. Windsor, The detection of buried pipes from time-of-flight radar data, IEEE Trans. Geosci. Remote Sens. 46 (8) (2008) 2254–2266. [18] C.G. Windsor, L. Capineri, P. Falorni, A data pair-labeled generalized hough transform for radar location of buried objects, IEEE Geosci. Remote Sens. Lett. 11 (2014) 124–127. [19] W. Li, X. Cui, L. Guo, J. Chen, X. Chen, X. Cao, Tree root automatic recognition in ground penetrating radar profiles based on randomized hough transform, IEEE Geosci. Remote Sens. Lett. 8 (5) (2016) 430. 1–16. [20] H. Zhou, M. Tian, X. Chen, Feature extraction and classification of echo signal of ground penetrating radar, Wuhan Univ. J. Nat. Sci. 10 (6) (2005) 1009–1012. [21] A. Dell’Acqua, A. Sarti, S. Tubaro, L. Zanzi, Detection of linear objects in GPR data, Signal Process. 84 (4) (2004) 785–799. [22] E. Passoli, F. Melgani, F. Donelli, Automatic analysis of GPR images: a patternrecognition approach, IEEE Trans. Geosci. Remote Sens. 47 (7) (2009) 2206– 2217. [23] X. Xie, H. Qin, C. Yu, L. Liu, An automatic recognition algorithm for GPR images of RC structure voids, J. Appl. Geophysics 99 (2013) 125–134. [24] H.S. Youn, C.C. Chen, Automatic GPR target detection and clutter reduction using neural network. in: 9th International Conference on Ground Penetrating Radar, S. Barbara, USA, 2002 [25] S. Shihab, W. Al-Nuaimy, Y. Huang, A. Eriksen, Neural network target identifier based on statistical features of GPR signals. in: 9th International Conference on Ground PenetratingRadar, S. Barbara, USA, 2002 [26] W. Al-Nuaimy, Y. Huang, M.T.C. Fang, V.T. Nguyen, A. Erikson, Automatic detection of buried utilities and solid objects with GPR using neural networks and pattern recognition, J. Appl. Geophysics 43 (2–4) (2000) 157–165. [27] A. Giannopoulos, Modelling ground penetrating radar by GprMax, Construct. Build. Mater. 19 (10) (2005) 755–762. [28] F. Frezza, L. Pajewski, C. Ponti, G. Schettini, N. Tedeschi, Cylindrical-Wave Approach for electromagnetic scattering by subsurface metallic targets in a lossy medium, J. Appl. Geophysics 97 (2013) 55–59. [29] S. Birkenfeld, Automatic detection of reflexion hyperbolas in GPR data with neural networks, World Automation Congress (WAC 2010), Kobe, Japan, 2010, pp. 1189–1194.

Please cite this article in press as: A. Ristic´ et al., Integration of modern remote sensing technologies for faster utility mapping and data extraction, Constr. Build. Mater. (2017), http://dx.doi.org/10.1016/j.conbuildmat.2017.07.030

16

A. Ristic´ et al. / Construction and Building Materials xxx (2017) xxx–xxx

[30] U.S. Khan, W. Al-Nuaimy, F.E.A. El-Samie, Detection of land- mines and underground utilities from acoustic and GPR images with a cepstral approach, J. Vis. Commun. Image Represent. 21 (7) (2010) 731–740. [31] C. Maas, J. Schmalzl, Using pattern recognition to automatically localize reflection hyperbolas in data from ground penetrating radar, Comput. Geosci. 58 (1) (2013) 116–125. [32] X. Núñez-Nieto, M. Solla, P. Gómez-Pérez, H. Lorenzo, GPR Signal characterization for automated landmine and UXO detection based on machine learning techniques, Remote Sens. 6 (10) (2014) 9729–9748. [33] M.R. Shaw, S.G. Millard, T.C.K. Molyneaux, M.J. Taylor, J.H. Bungey, Location of steel reinforcement in concrete using ground penetrating radar and neural networks, NDT&E International 38 (3) (2005) 203–212. [34] P. Gamba, S. Lossani, Neural detection of pipe signatures in ground penetrating radar images, IEEE Trans. Geosci. Remote Sens. 38 (2) (2000) 790–797. [35] A. Plan, Ground Penetrating Radar- and Resistivity Measurements – a Surveys of Corrosion Enhancing Soil Properties Around District Heating Pipes in Ängelholm Dissertations in Geology at, Lunds university, 2016. 34 p. [36] C. Persson, Predicting of Long-Term Insulation Performance of District Heating Pipes, Chalmers University of Technology, Göteborg, Sweden, 2015, p. 73p. [37] O. Friman, P. Follo, J. Ahlberg, S. Sjökvist, Methods for large-scale monitoring of district heating systems using airborne thermography, IEEE Trans. Geosci. Remote Sens. 52 (8) (2014) 5175–5182. [38] T. Ahonen, A. Hadid, M. Pietikainen, Face Recognition with Local Binary Patterns, University of Oulu, Finland, 2004. [39] N. Dalal, B. Triggs, Histograms of Oriented Gradients for Human Detection, IEEE Computer Society Conference on Computer Vision and Pattern Recognition, France, 2005 [40] W. Al-Nuaimy, Automatic feature detection and Interpretation in groundpenetrating Radar data, University of Liverpool, UK, 1999.

Aleksandar Ristic received the PhD degree in Electrical engineering and geo-informatics in 2009 at Faculty of technical sciences, University of Novi Sad. He is currently an Associate Professor in the Department of automation, geomatics and control systems and lecturer in underground utility detection by GPR and EML, geosensor networks and control systems in geomatics. His current research interests include development of innovative inspection procedures for GPR surveying of underground utilities in urban areas, and quantitative estimation of EM and physical properties from GPR data with development of advanced GPR data processing techniques. Dr. Ristic participated in the project of development and implementation of first GPS permanent stations network in Serbia, in 2003. He published a number of papers in journals and scientific conferences proceedings related to GPR technology. Also, he is reviewer in several journals. Dr. Ristic participated in international projects (EUPOS Berlin, COST 1208). He has been project leader in a number of projects directly related to GPR applications and runs the laboratory for sub-terrestrial remote sensing.

ˇ eljko Bugarinovic´ Received his M.Sc. degree in 2014 Z from the Faculty of technical sciences, Novi Sad, study program Geodesy and Geomatics. Currently is Ph.D. student of the Faculty of Technical Sciences in Novi Sad. He works as assistant at Department of automation, geomatics and control systems and teaches exercises in subjects related to application of GPR, GNSS and other geoinformation technologies. His area of research are non-destructive technologies (GPR and other) for underground utility detection, algorithms for automatic data extraction, pattern recognition and remote sensing.

Milan Vrtunski received the M.Sc. degree in Electrical engineering and computer science in 2006 at Faculty of technical sciences, University of Novi Sad. He is currently Ph.D. student. He works as assistant at Department of automation, geomatics and control systems and teaches exercises in subjects related to application of GPR, GNSS and other geoinformation technologies. His current research interests include development of advanced methods for joint data acquisition by GPR and other ND technologies, data processing and application of GPR for determination of soil structure with volumetric moisture content. Mr. Vrtunski published a number of papers in journals and scientific conferences proceedings and participated in many projects directly related to GPR applications.

Miro Govedarica received the PhD degree in Geoinformatics and information systems in 2001 at Faculty of technical sciences, University of Novi Sad. He is currently a Full Professor in the Department of automation, geomatics and control systems and lecturer in geoinformatics, remote sensing and geospatial services. His current research interests include development of innovative inspection procedures for GPR surveying of underground utilities in urban areas, development of advanced GPR data processing techniques and development of advanced methods for joint processing of data collected by GPR and other non-destructive testing techniques. Dr. Govedarica is founder (2003) of the Geospatial Technologies and Systems Center (GTSC). He participated in the project of development and implementation of first GPS permanent stations network in Serbia, in 2003. He published a number of papers in journals and scientific conferences proceedings. Also, he is reviewer in several journals. He participated in number of projects directly related to GPR applications, and in international projects (EUPOS Berlin, COST 1208 MC Member) as well. Dr. Govedarica is Head of the studying program Geodesy and geomatics at Faculty of technical sciences.

Dusan Petrovacki received his PhD degree in 1979, at Faculty of technical sciences, University of Novi Sad. He is IEEE society member. His current research interests include applications of geospatial and information technologies in various fields, quantitative estimation of EM and physical properties from GPR data with development of advanced GPR data processing techniques. Dr. Petrovacˇki is Head and founder (2003) of the Geospatial Technologies and Systems Center (GTSC) and former Dean of Faculty of Technical sciences (19911998). Also, he was Regional secretary (Vojvodina/Serbia) for Science and Technology (2002-2004)., Novi Sad. He lead the project of development and implementation of first GPS permanent stations network in Serbia, in 2003. He participated in 10 international projects (UNDP, MSF (US), DOE (US), Academia Sinica (China), British Council, Erasmus University Rotterdam (Netherlands), EUPOS Berlin). He published more than 150 papers in journals and scientific conferences proceedings.

Please cite this article in press as: A. Ristic´ et al., Integration of modern remote sensing technologies for faster utility mapping and data extraction, Constr. Build. Mater. (2017), http://dx.doi.org/10.1016/j.conbuildmat.2017.07.030