The extensive demand for the rapid and high quality acquisition of 3D urban features has motivated much research effort towards increased automation. Road network database development is expensive and time consuming since most of the road network extraction depends heavily on human labor. Creating, maintaining, and updating the transportation network database is an essential requirement for many applications such as automated vehicle navigation, traffic management, safety, disaster assessment and future planning. Moreover, many critical decisions such as emergency response, evacuation, and incident management could be based on such data bases (Xiong, 2001).

With the availability of many sources of data such as conventional imagery, Synthetic Aperture Radar (SAR) imaging, Interferometric Synthetic Aperture Radar Digital Elevation Models (IFSAR DEMs), and LiDAR DEMs, there are many avenues open to derive terrain and feature data in urban areas. Through much research, it has been shown that laser scanning data has the potential to support 3D feature extraction, especially if combined with other types of data such as 2D GIS ground plans (Maas, 1999; Brenner and Haala, 1999; Weidner and Förstner, 1995). Despite the fact that LiDAR data is attractive in terms of cost per high quality data point, the quantity of the data makes a challenge for storage and display (Maas and Vosselman, 1999). Acquiring 3D object descriptions from such data is a challenging problem and many approaches have been tried to solve it. Several of them have succeeded with some limitations. The principal idea of this research is to extract the road network from laser data exclusively, utilizing both laser intensity returns and height information.

In much previous work, road extraction research activities were more concentrated on extracting roads from imagery in rural areas than urban areas due to the complexity of urban scenes (Hinz and Baumgartner, 2000). Moreover, road extraction from aerial imagery has some limitations due to several factors. Reflectance dependency on viewing conditions and illumination, variability of the surface material, occlusion and resolution limitations are some examples of those factors (Price, 1999).

Recently, pattern recognition and a hierarchical road model have been utilized to improve the extraction procedures and consequently the results. A detailed and scaled road model also was used to update road maps (Vosselman and Gunst, 1997). This method works very well in open rural areas but it is sensitive to shadows and occlusions.

The road detection procedure described here includes first segmenting roads from other natural and man-made features such as trees and buildings. The aim is to develop a fast and accurate method to detect the road network in urban areas from LiDAR data. Filtering to identify candidate roads from other urban features was the first step in this work. Many segmentation techniques are feasible with LiDAR data such as direct thresholding determined by histogram analysis, the use of 2D GIS data, and multispectral and hyperspectral inferences (Maas, 1999; Brunn and Weidner, 1997). Direct thresholding works by extracting the terrain surface using a filtering technique such as the morphological filter. The segmentation of roads from other urban features in this research was based on the resulting intensity from the laser scanning survey. However, due to some similarity in the reflections between some features, LiDAR heights were also necessary to support the segmentation procedure and to refine it.

In this work, the intensity was used to segment roads from other objects in the data set. However, as in many image interpretation tasks, the result of segmentation is not always complete and this problem is more obvious in road extraction in particular. So, more steps were needed to refine this classification. Besides the intensity analysis, the LiDAR heights and their inferences was used as a key to discriminate between ground and non-ground features since roads are known to be ground features with few exceptions such as in bridges. The DTM (Digital Terrain Model) was extracted by applying a modified mathematical morphological filter on the last return heights of the digital surface model (LDSM). Consequently, non-ground features were detected using the normalized LDSM (the difference between LDSM and DTM), and replaced by the corresponding terrain. Then, this result was utilized to filter out misclassified pixels from the previous step. This approach was tested on the data set that has been collected over the Purdue University campus in spring 2001 with an approximate density of one data point per each square meter.

Data Classification Procedure

Laser point clouds must be segmented in order to distinguish between building, roads, trees and other urban features. It is a mandatory process to differentiate among diverse objects in the scene. Extraneous objects such as trees that do not belong to the road category should be detected and removed from the scene. The filtering process, which was used in this project to segment roads from other objects, consists of two main steps. The first one is to utilize the return intensity image of the LiDAR data in order to detect roads based on their surface reflectivity. But due to similarity in feature reflectance another step was needed to filter out the undesired features. The second step is making use of "first minus last" return analysis and utilizing the 3D LiDAR information inferences to filter out misclassified pixels that do not belong to road category.

Intensity Return Analysis

Most laser scanning systems today have the ability to record the intensity of the reflected laser pulse. The data that we have tested was collected using an Optech ALTM1210 LiDAR system, which was operated by Woolpert Consultants, and it has the ability to record the laser reflection intensity. This particular system operates at wavelength of 1.064 µm, i.e., in the near infrared region of the spectrum. In current ranging systems, the laser wavelength can occupy any segment of the optical band between 0.8µm and 1.6µm. However, the selection of the optical wavelength of the laser depends on the design aspects of the laser system and the purpose of the survey. The way to measure the intensity of a reflected pulse is by quantizing the amplitude when receiving the echo. This technique is used in most current laser sensors. It produces a noisy image since it does not reflect the intensity of the whole span of the received echo, rather only an instantaneous intensity. However, the intensity image can be used to improve the laser scanning data classification and filtering. Beside other factors such as wavelength and flight height, the reflectance varies from one scanned surface to another based on the surface properties. An object surface like sand has a reflectance range between 10% and 30%, vegetation between 30% and 50%, while snow and ice can have up to 80% reflectance (Wever and Lindenberger, 1999). Calm water surfaces introduce strong dependence on view angle since they become specular rather than diffuse surfaces. Dark surfaces such as asphalt or slate roofs absorb the laser beam and consequently no reflection will be recorded which results in a gap in LiDAR data.

The principle of this approach is to start with the segmentation of the intensity image of the LiDAR data. First, the histogram of the digital numbers from the quantized intensity over the tested area was analyzed. It was found that the quantized intensity values range from zero to over 8000 with a mean of 330 and a standard deviation of 230. However, more than 95% of the values range from zero to 800 as shown in Figure 1. Those response values are based on a scale starting at zero which represents the weakest return in the data set. By analyzing the data we found that high intensity (larger than 800) are small in number and they are mostly positioned in parking lots. After inspecting the intensity image closely it has been found that those high intensity spots belong to specular effects on glass and metal parts of cars, and occasionally, buildings. Those high responses were then filtered out based on the histogram interpretations of the desired range and replaced by the maximum acceptable range, which is 800 in our case. The histogram of the filtered intensity responses is shown Figure 2.

Since the raw LiDAR data is in an irregular form, the intensity returns were interpolated and resampled into a regular structure in order to view and process the intensity image. The resampling was done with a spacing of one meter since the irregular spacing was approximately one meter. Figure 3 shows the filtered and scaled (0-250) intensity image. The following step was to define the intensity range that road pavement materials might have in the data set in order to supervise the classification procedure. To accomplish that, training samples were chosen based on prior knowledge of road locations to estimate the reflectivity range of roads in the tested area. The training samples were chosen carefully from different locations among the data set (as shown in Figure 3) to cover the whole range of pavement materials. Consequently, the decision boundaries were obtained for the road class. It was obvious from the intensity image (Figure 3) that road pavement materials have low reflectance. According to the training samples the intensity values of roads range between 18 and 60. Accordingly, the intensity image was classified into two classes, roads and non-roads, based on the intensity returns for each resampled pixel. The classification result is shown in Figure 4.

However, due to the similarity in reflectance between objects such as asphalt road and building roofs, which maybe be coated with similar materials, there was a significant misclassification error in this procedure. Also sidewalks and some areas of vegetation have similar returns to roads. As a result of that, a refinement step was needed. This step utilizes the LiDAR 3D information, which is the main part of the LiDAR data, as discussed in the following section.

Classification Refinement

The classification refinement step was designed to address the misclassified areas and improve the classification accuracy. The misclassification was mostly a positive error; i.e., areas were classified as roads while they are not in reality. Most of the misclassified areas are buildings roofs, sidewalks, driveways, and grass. Therefore it was mandatory to detect those objects in order to improve the classification result. In order to detect grass and other ground vegetation, this step utilizes the capability of the LiDAR system of being able to measure and record two returns (first and last) from each pulse (Alharthy and Bethel, 2002). Also the interpretations of laser heights and extracted DTM provide good indication of the position of an object, i.e., whether it is a terrain or above terrain object. Those two algorithms and the way they were utilized are discussed in brief below.

First-Last return Analysis

The LiDAR system that was used to collect this data was an Optech ALTM1210 and it has the ability to capture two returns (first and last) per each height point. This is due to the fact that the laser pulse is not a single ray but occupies an extended solid angle. It has an angular beamwidth and its footprint will a take circular shape when it reaches the ground. Based on both laser and scene characteristics and physics, the laser beam could penetrate some objects. Therefore some of its energy will be reflected back from the object top surface and others might penetrate to different depths before they are reflected. Generally, this produces an extended return signal. The first and last will be recorded and their corresponding heights will be calculated separately. As a result of that, two heights were assigned to each data point. The first return heights contain more noise since they will reflect every object on the ground such as trees, bushes, cars, etc. On the other hand, the computed heights based on the last received return represent only the non-penetrable objects such as the ground, and buildings. Two different heights of one point give an indication of the presence of a penetrable object such as trees and vegetation or possibly "straddling" a multi-height object. On the other hand, if a data point has the same height for first and last returns, then this point belongs to a non-penetrable 'solid' object. This step is just to locate the penetrable objects, which are categorically considered non-road regions.

A high level of discrepancy between first and last return heights takes place only on penetrable objects, which is considered as noise, as stated above, according to the objective of this work. The discrepancy map was then employed to improve the classification. If a pixel among those that were classified as road from the first step has a discrepancy larger than a certain threshold, this pixel will be considered as a penetrable object and consequently a non-road pixel. Accordingly, it will be filtered out. The result of this refinement step is shown in Figure 5. We can see how the trees on roadsides and grass in open areas are filtered out. During this step the classification was improved as shown in Table 1. However, there are some parts of the misclassified areas still not filtered out. Those parts are mostly buildings roofs that are coated with materials similar to road pavements materials.

3D Information and Terrain Non-terrain Point Analysis

As stated above, some of the roof buildings were covered with similar materials to those that were used in road pavement. This similarity causes both surfaces to have comparable intensity returns when they are illuminated with the laser pulse. Consequently such roofs were classified as roads during the intensity segmentation procedure. Moreover, they are not filtered out in the first-last return analysis since roofs are non-penetrable objects. Accordingly, in order to reclassify those roofs "non-ground objects," they need to be detected first. The direct route to distinguish between terrain and above terrain objects is to extract the terrain model first. Lindenberger (1993) introduces and shows how mathematical morphology can be used for filtering data recorded by a laser profiler and consequently for DTM generation. (Vosselman 2000) introduces an approach called slope based filtering in which he established a modified version of the mathematical morphology. In this work, a similar mathematical morphology filter in an autoregressive mode was used to extract the DTM as shown in equation 1. I represents a pixel in the given intensity image, G represents the computed gradient at that pixel using the last return heights, DSM is the raw resampled digital surface model, DTM is the extracted digital terrain model which represents the bare ground, and tG is the maximum allowed gradient and can be computed based on surface nature and scene characteristics. Then, each classified road pixel will be tested to determine weather it belongs to the bare ground model (DTM) or not by comparing its last LiDAR height with the extracted DTM. Accordingly, if it is not a ground pixel then it will be reclassified as a non-road one. The result of this refinement step was significant as shown in Figure 6 and Table 1.

DTM = {I (i,j) ? DSM ? ;I (i,j) :G(i,j) < tG}

Connected Component and Final Result

Connectivity between pixels is a useful concept used in obtaining homogenous regions in an image. So, to decide that two pixels belong to one region they should be adjacent and homogenous (have similar gray level for example). This concept was utilized here to establish road regions and consequently construct the road network. Also, the size of those regions was used to filter out small regions like side walks and other scattered noise. The result of this step is shown in Figure 7. Table 1 shows the summary of this procedure's performance.

Discussion

The aim in this work is to design a simple and fast method to detect road networks in urban areas using LiDAR data exclusively. We restricted the procedure to not require any other source of data other than LiDAR. This was done intentionally to avoid the limitation of availability of other sources of information in some areas. Sources such as ground plans, imagery and multispectral data are not available for every desired site and might be not up to date. This method starts with filtering the raw LiDAR data to remove "noise" unrelated to the road class. Refinement steps followed to improve the heuristic classification of the data. The performance of those steps is good, however, accuracy of the extracted road network needs to be investigated more. The presented algorithm works very well with the tested data and is expected to do the same with others. Yet, it might fail in some other areas where road materials vary within small areas and consequently their reflection will cover a wide range of the intensity scale. This surely will affect the classification procedure and increase the misclassification error. Practically, in general, road detection algorithms are able to work only with a limited set of characteristics. And those characteristics can be modified based on the data nature and specification. When these characteristics change beyond expected limits the algorithm may fail.

The quality of the result depends on both the algorithm and the data. In general, LiDAR data is not continuous and it does not represent the entire scene since there are some gaps between pulse footprints on the ground. Thus we need to treat the data as sampling of a continuous surface with an aperture corresponding to the ground footprint of the laser. With a sampling interval of one meter, ideally to avoid aliasing, the aperture would also be one meter. In this case, the aperture (laser footprint) was about 30cm so the data was undersampled to that degree.

Therefore, as shown in Figure 7, the resulting road network is not totally complete due to some classification errors as stated above. However, these results seem to be promising compared with other data sources and approaches. These results are preliminary and more development of the algorithm needs to be performed.

Ackermann, F. (1999). Airborne laser scanning-present status and future expectations. ISPRS Journal of Photogrammetry & Remote Sensing, 54:64-67.

Alharthy, A., Bethel, J. (2002). Heuristic filtering and 3d feature extraction from LiDAR data. ISPRS Commission III, September 9-13, 2002, Graz, Austria.

Brenner, C. (2000). Towards fully automated generation of city models. ISPRS, vol. XXXIII, Amsterdam 2000.

Brenner, C., Haala, N. (1999). Extracting of building and trees in urban environments. ISPRS Journal of photogrammetry & Remote Sensing, 54:130-137.

Brenner, C., Haala, N. (1999). Rapid production of virtual reality city models. Geoinformationssysteme, (2):22-28.

Brunn, A. (2001). Statistical interpretation of dem and image data for building extraction. In A. Grün, editor, ASCONA-Workshop 2001. Balkema-Verlag.

Brunn, A., Weidner, U. (1997). Extracting buildings from digital surface models. IAPRS Vol. 32, part 3-4w2, Stuttgart.

Hinz, S., Baumgartner, A. (2000). Road Extraction in Urban Areas supported by Context Objects. In: International Archives of Photogrammetry and Remote Sensing, Vol. 33, part B3.

Hug, Ch., Wehr, A. (1997). Detecting and identifying topographic objects in imaging laser altimetry data. IAPRS 32, 19-26, part 3-2w3.

Jelalian, A. (1992). Laser Radar System, Artech House, Boston.

Lindenberger, J., (1993). Laser-Profilmessungen zur topographischen Geländeaufnahme. Deutsche Geodätische Kommission, Series C, No. 400, Munich.

Maas, G-H. (1999). The potential of height texture measures for the segmentation of airborne laserscanner data. Presented at the 4th Airborne Remote Sensing Conference and Exhibition, Ottawa, Ontario, Canada, 21-24 June 1999.

Maas, H., Vosselman, G. (1999). Two algorithms for extracting building model from raw laser altimetry data. ISPRS Journal of Photogrammetry and Remote Sensing, 54(2-3):153-163.

Mikhail, Edward M., Ackermann, F. (1976). Observation and Least Squares.

Price, K., (1999). Road grid extraction and verification. In: International Archives of Photogrammetry and Remote Sensing, Vol. 32, part 3-2W5, pp. 101-106.

Song, J., Han, S., Yu, K., Kim, Y., (2002). Assessing the Possibility of Land-Cover Classification Using LiDAR Intensity Data. ISPRS Commission III, September 9-13, 2002, Graz, Austria.

Uwe Weidner, Wolfgang Förstner (1995). Towards Automatic Building Reconstruction from High Resolution Digital Elevation Models. ISPRS-Journal, 50(4), pp.38-49, 1995.

Vosselman, G., Gunst, M., (1997). Updating road maps by contextual reasoning. Proceedings 2nd Workshop on Automatic Extraction of Man-Made Objects from Aerial and Space Images, Ascona, May 5-9.

Vosselman, G., (2000). Slope based filtering of laser altimetry data. IAPRS, Vol XXIII, Part B3, pages 935-942. Amsterdam, The Netherlands.

Xiong, Demin. (2001). Automated road network extraction from high resolution images. Technical notes, National consortium for safety, issue 3, may 2001.