INTEGRATION OF LIDAR DATA INTO A MUNICIPAL GIS TO STUDY SOLAR RADIATION

Identifying the right roofs to install solar panels inside a urban area is crucial for both private citizens and the whole local population. The aim is not easy because a lot of consideration must be made: insolation, orientation of the surface, size of the surface, shading due to topography, shading due to taller buildings next the surface, shading due to taller vegetation and other possible problems typical of urban areas like the presence of chimneys. Accuracy of data related to the analyzed surfaces is indeed fundamental, and also the detail of geometric models used to represent buildings and their roofs. The complexity that these roofs can reach is elevated. This work uses LiDAR data to obtain, with a semi-automatic technique, the full geometry of each roof part complementing the pre-existing building data in the municipal cartography. With this data is possible to evaluate the placement of solar panels on roofs of a whole city analyzing the solar potential of each building in detail. Other traditional techniques, like photogrammetry, need strong manual editing effort in order to identify slopes and insert vector on surfaces at the right height. Regarding LiDAR data, in order to perform accurate modelling, it is necessary to obtain an high density point cloud. The method proposed can also be used as a fast and linear workflow process for an area where LiDAR data are available and a municipal cartography already exist: LiDAR data can be furthermore successfully used to cross-check errors in pre-existent digital cartography that can remain otherwise hidden.


INTRODUCTION
A deep knowledge of the city and urban structure geometry is today essential for a great number of applications.Buildings shape and their placement are fundamental informations needed to conduct a better management of the city itself.Today modern geomatic offers efficient tools with high precision to obtain data on a larger scale.Three-dimensional representation in form of 3D city models can be realized with different techniques, in particular digital photogrammetry and laser scanning that can output products at different level of details and precision, with great potential and the possibility to use into Geographic Information Systems maintained by municipalities and public agencies for urban management and planning.Knowing that a city is a central hub of energy use, among useful applications there is the study about energy consumption and production.These analysis can be conducted from a urban scale even to single buildings in detail.One example is Bologna Solar City, from Municipality of Bologna (Italy): a web application created to evaluate the placement of solar panels on roofs.It will be discussed in detail in the following section.Another example is to use LiDAR data and a derivative Digital Surface Model to produce a true ortophoto (Günay et al., 2007) that is useful in order to refine a detailed energy map relative to each building.This kind of analysis at building level was performed within EnergyCity, a project from Central Europe, to better understand the actual energy management and support strategies from the urban municipalities in order to reduce energy consumption and CO2 emissions (Bitelli and Conte, 2011).

Bologna Solar City
Municipality of Bologna has a long experience in cartography and a good collection of data over the years in digital format.These data are continuously updated through different sources: georeferenced CAD drawings from executive projects, manual extraction of features from aerial orthophoto taken every year and surveys on the field are used to update the original large scale numerical cartography (CTC) made by aerial photogrammetry.In 2010, in order to help its inhabitants to reduce energy costs and save the environment, Bologna SIT office created a web application called Bologna Solar City to search for renewable energy systems applied in the Bologna territory and to evaluate the placement of solar panels on roofs.The basic idea was to estimate the solar energy available on roofs.This evaluation is the sum of direct and diffuse solar energies and it was calculated with Ar-cGIS tools.The algorithm applied required several elements in input, first of all a DSM (Digital Surface Model) of the city territory developed with many altimetric information deducted from the municipal cartography.The other factors considered are, on a raster-based model, weather conditions, position and shade of each location (Minghetti et al., 2011).The map services were published using ArcGIS Server and a user friendly interface application was made available using ArcGIS API for Flex.Through this web application it is possible to calculate the total solar radiation per building, thus estimating the potential yearly energy production and the reduction of CO2 emissions.
International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XL-1/W1, ISPRS Hannover Workshop 2013, 21 -24 May 2013, Hannover, Germany The main simplification in Bologna Solar City was to consider all the roofs as flat, because in 1:2000 municipal cartography there are only two height information for each building: the bottom line and the height of the rain gutter.Unfortunately the greatest part of roofs in Bologna, as well as in many other cities, are composed from a complex arrangement of different slopes, with many ridges and valleys.In order to perform an accurate analysis, different parameters must be considered: insolation, orientation of the surface, size of the surface, shading due to topography, shading due to taller buildings next to the surface, shading due to taller vegetation and other possible problems typical of urban area like the presence of chimneys.
Accuracy of data relative to analyzed surfaces is indeed fundamental, and also the detail of geometric models used to represent buildings and their roofs.For this reason the precise geometry has to be acquired in order to have a good evaluation about the solar potential of each part of the roof.To obtain this kind of information a Digital Surface Model of the city could be very useful.This is obtained, in the presented work, from LiDAR data using point clouds acquired from airborne laser scanning.These points are loaded in ENVI LiDAR, a software from Exelis Visual Information Solutions.This software can classify and also directly extract 3D features from the point cloud, if the data density is adequate.One challenge is to integrate LiDAR data and digital 2D cartography taking in account that each dataset has a different origin, a different reference system and also a different acquisition epoch.

AVAILABLE DATASETS
At the moment of writing this paper there is no LiDAR data available to cover the whole city of Bologna, but only partial coverage.
The chosen study area is Reno district, one of the areas with full data coverage along with Porto district which is less interesting in the analysis of this paper, being partially covered by railway and Station of Bologna and being partially inside the historical center of the city, where the installation of solar panels on roofs has to follow strict rules.The dataset received is in LAS format in a projected CRS: WGS84 / UTM zone 32N (EPSG:32632).Elevation data is therefore referred to WGS84 ellipsoid.
The LiDAR cloud is not particulary dense.The survey was not performed with 3D building extraction as primary goal.Therefore, the average data density is only equal to 1 point/m 2 with higher value where two adjacent strips cover the same area as can be seen in Figure 2.
In this work were also used the shapefiles from digital GIS cartography (CTC) at a nominal 1:2000 scale, produced and mantained by Municipality of Bologna and now shared as Open Data: • Buildings perimeter The result is more accurate than transformation parameters inside GIS software.

DATA PROCESSING
The study here described uses the software ENVI LiDAR from Exelis VIS.The procedure firstly consists in importing all the 20 LAS files acquired over Bologna in ENVI LiDAR that converts them into a unique binary data file for a faster manipulation.The first operation to perform is to check data density through the work area.This is illustrated in Figure 2 exporting the result as a raster.
ENVI LiDAR procedures operate to filter data and classify each point of the cloud.Then, through automated feature identification with manual tools for quality control, it is possible to extract different features as: vegetation, power line vectors and power poles, DSM and filtered DTM, building vectors with roof face and perimeter.
One relevant improvement over traditional LiDAR software is in further analysis after point classification.With appropriate segmentation procedures it's possible to create vectors and then interpolating surfaces among planes defined in the point cloud (Rottensteiner and Briese, 2002).The software, with a semiautomatic procedure among building's vectors, can identify the correct position, aspect and slope of each roof plane in the work area extracting consistent 3D building models geometrically correct.

Fine-tuning parameters
Regarding building extraction some production parameters has to be set: • Minimum area -building planes with an area below set value, will not be vectorized.
• Near Ground Filter Width -filters out buses, trucks and train cars.Building planes less than 5 meters from the ground are not classified as buildings unless they have a width greater than this parameter.
• Buildings Points Range -spatial variation of building's points: the value is used by the software's algorithm that scans for planar surfaces when the point density is not constant inside analyzed area or when there are some holes in the point cloud dataset.In this cases the value can be set as high as 1.2 or 1.4 meters, increasing risk to misclassify vegetation point into building class.
• Plane Surface Tolerance -curved roofs are modeled using a series of successive planes.A new plane starts when the distance between the points being analysed and the previous plane reaches this value.Plane Surface Tolerance (PST) is a sensible parameter.The default value is 30 cm.With the available data best results can be achieved using a value of 45 cm as seen in Figure 3, due to point cloud low density and low accuracy.However, using such an high value has a disadvantage when it comes to modelling curved surfaces as can be seen in Figure 4a, where a curved roof of an industrial building is extracted.
Within the urban environment where the study area is set, there are various problems related to features, like buildings from various periods with hipped or gabled roof often irregular and made more complex by the presence of heterogeneous elements.There is also a strong noise effect on point cloud data caused by vegetation that is located adjacent to the buildings on the fringes of urban roads.In correcting the slopes of the roofs the primary focus must be to keep their correct aspect rather than the exact representation of all the complexity that some buildings show as chimneys, dormers, and complex geometries.It would not be technically possible to represent them due to the reduced density of the points, it is also not relevant to the objectives that have been placed in this research, which require to take greater account regarding orientation of the surfaces analyzed with the aim of evaluating their solar radiation potential.
Through extensive tests the following parameters were chosen to perform analysis on the dataset for the test area: • Minimum area: 5 m 2 • Near Ground Filter Width: 150 cm • Buildings Points Range: 1,4 m • Plane Surface Tolerance: 45 cm An example of the results is shown in Figure 5.

Elevation data
With airborne laser scanner a good accuracy on vertical position can be obtained, useful to create terrain models, but the accuracy on planimetric position is generally worse and this is relevant for analysis like mapping of buildings as proposed in this work (Maas, 2003).Regarding planimetric accuracy of LiDAR data there is a main influence by the scanning mirror mechanism and the laser spot diameter which is directly proportional to altitude Above Ground Level.To obtain a good accuracy is better to mantain low flight height and set laser beam divergence to minimum in order to have the smallest spot possible at ground level.Available data has a low horizontal accuracy, which in this case is 0.6 meters, about 1/2000 flight height.With modern instruments or flying at lower altitude a better result can be of course achieved.
International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XL-1/W1, ISPRS Hannover Workshop 2013, 21 -24 May 2013, Hannover, Germany Regarding height accuracy, tests were performed with ConvER3 GPS7 and data from the official IGM grid in the whole Reno district.The Geoid Height values in Reno district are between the range from 39.2 meters to 39.4 meters.The height from Elevation points in shapefile from CTC is orthometric while the data at the same coordinate in LiDAR cloud is referred to WGS84 ellipsoid.The local Geoid Height, computed as a difference between the two values, is close to the corrispondent IGM value considering the accuracy of both LiDAR and shapefile data.
It's also possible to use LiDAR data to cross-check possible errors in pre-existent digital cartography, that can remain otherwise hidden, as can be seen in Figure 6.One kind of analysis can be to highlight areas where the difference between known Geoid Height value in the area and the calculated Geoid Height (difference between ellipsoid height from LiDAR and orthometric height from numerical cartography) isn't close to the given tolerance.
Figure 6: Possible error highlighted with a yellow triangle, Orthometric Height from CTC in green, calculated Geoid Height in red

REFINEMENTS, EDITING AND QUALITY ASSURANCE
Using ENVI LiDAR after the beginning classification procedure and the further building extraction procedure, is then possible to use tools to ensure that resulting feature representations are accurate.In the Quality Assurance mode is possible to interactively fine-tune point classifications and specially modify feature vectors such as building's roof shape using various tools such as cross sections.Some points do not belong to building class, as can be seen in Figure 7. Laser scanners detects points on antennas, chimneys and every reflective surface.These points and other outliers are automatically filtered out from software algorithm as they become noise around the roof surface.With tools like Cross Section in ENVI LiDAR is possible to manually remove these points if they were not correctly classified.Also, with Change Classification Code tool, points that were not classified as buildings can be assigned to the correct class.If the data density is adequate, it is also possible to set the selected points as a starting point for a new interpolating surface searching.This reduces limitations set on noise in the automatic process, thereby increasing the ability of ENVI LiDAR to set a plane at the designated point.Each vector can be manually edited, moving its vertices in every direction, to fit better to the local point cloud (Figure 8).A red triangle pointing down will appear highlighting the area where the modification is applied.All the possible operations to be performed on the data are listed in a dialog box.

Figure 8: Manual editing procedure
There is also the chance to introduce external data into the point cloud.This is a useful feature in order to verify the correct placement of buildings and if some of them have initially remained unnoticed from the automatic 3D extraction.The municipal numerical cartography is loaded to perform a visual check.As can be seen in Figure 9 some buildings were not recognized from the low density point data cloud.To perform a numerical quality-check against another data source, the extracted footprints from LiDAR data can be compared with the footprints from pre-existing digital cartography using a confusion matrix and Kappa Index (Shirowzhan and Lim, 2012).In order to perform this comparison, it must be considered that in many cases LiDAR data can identify roofs areas, but especially where there is a pitched roof that extends over perimetral walls, the extension of that surface projected at ground level is wider than the area delimited by building's walls.
International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XL-1/W1, ISPRS Hannover Workshop 2013, 21 -24 May 2013, Hannover, Germany Furthermore, in a traditional survey or in photogrammetry it is possible to identify quality-control elements in corrispondence with surface discontinuity as walls intersection or rain gutter line in roofs.Those two kind of controls are not available in LiDAR data especially when there is a low density in points.In pitched roofs, the clearest discontinuity is the ridge roof line which position is determined at the intersection of the adjacent roof planes.
An evaluation of planimetric accuracy can be performed comparing ridge roof lines in overlapping LiDAR strips (Vosselman, 2008).Those lines can be easily extracted from the data cloud with software like ENVI LiDAR.
The 3D models of buildings superimposed on the photogrammetric model in a stereophotogrammetric workstation can provide the way for a quick visual inspection of the correct description of buildings geometry.Once the interested buildings passed Quality Assurance control (Figure 10) the polygon file is ready to export.As it's shown in Figure 11, the final 3D model can be viewed in a realistic representation, allowing to further qualify the accuracy of results flying through the scene.In order to perform solar radiation analysis, one possibility is to directly use the Digital Surface Model produced from LiDAR data and buildings shapefile from digital cartography as mask for the DSM to analyze only roofs area.The LiDAR points must be first classified in order to remove vegetation points and possible outliers.Doing so, the DSM can be computed only on points belonging to buildings and ground level.But with data so low in density as 1 point/m 2 the resulting DSM has a potential cell size of only 1 meter.Performing slope and aspect analysis at building level with so low detail brings in a lot of noise and salt and pepper effect that makes the data not suitable to identify every different surface of each roof.As can be seen in Figure 12 there are also holes in the data that occur where the point cloud density goes below 1 point/m 2 .
Figure 12: Aspect Map generated from DSM masked with Municipality shapefile To have a better resolution the point density must grow exponentially: the operative rule to obtain a P pixel size in generated DSM is P=1/sqrt(d) where d is the point cloud density.Using a dataset with 4 points/m 2 the raster resolution can be 0,50 m; using a dataset with 16 points/m 2 a raster resolution of 0,25 m can be achieved (Beinat and Sepic, 2008).
To obtain slope and aspect values in ESRI ArcGIS environment, using a 3D shapefile (Figure 13a) exported as a Polygon Z(M) from ENVI LiDAR, the best option is to transform 3D data to a raster.Generating a DSM allows to run 3D analysis tools like: line-of-sight, viewsheds or solar radiation calculations.This transformation is performed creating first a Triangular Irregular Network (TIN) using the roof polygons as breaklines and using constrained Delaunay triangulation where no densication occurs and each breakline segment is added as a single edge (Figure 13b).Then it's used a tool to generate a raster linearly interpolating cell z-values from the input TIN at the specified resolution of 0.5 meters that should be adequate given the instrument accuracy (Figure 13c).Then Extract by Mask function is performed to delete values outside the modeled buildings, using the building's perimeter polygons obtained from ENVI LiDAR as a mask (Figure 13d).With this final raster is possible to create derivative products like aspect (Figure 13e) and slope maps and use the resulting DSM to perform solar radiation analysis.
The ArcGIS Area Solar Radiation tool computes the global solar radiation on a surface for any time of day for each month of the year, investigating in every direction whether there are any obstacles, natural and artificial, and cross referencing this information with the map of the sun and the sky map.The first is based on the position of the sun and calculates direct beam radiation; while the second is based on the atmosphere and calculates diffuse solar radiation.Reflected solar radiation is excluded from the calculation (Minghetti et al., 2011).In order to perform a comparison between the two different solar radiation values, sloped roof against flat roof, the tool Set as Boxed was used: it calculates the perimeter of the roof at the rain gutter level and sets the roof to be flat as can be seen in Figure 14b.The global value obtained for the sample flat roof was compared with the solar radiation value among  As expected, the flat roof simplification tends to overstimate the solar potential of the building, while with proposed method that value can be computed with a good accuracy for each building.

CONCLUSIONS
The availability of software packages that allow to extract, virtually in automatic mode, 3D features, and in particular buildings, from high density LiDAR data is very attractive for a variety of applications in urban areas, among which regarding assessments in the energy field.A precise 3D roof's model can be easily extracted and with it all the information about surface area and orientation of each slope.
The work has shown that a software of this kind can be in particular useful in studies devoted to evaluation of solar systems installation on the rooftops (the example developed at Bologna Municipality has been described), or in analyses related to energy loss.Another secondary positive result is that this approach can help further validation of numerical cartography, derived from photogrammetric database, adopted as the backbone of a municipal GIS.
One of the main problems emerged is that, in order to reduce to a minimum the manual intervention in the presence of complex geometries of the roofs, it is important to have dataset with a density much greater than that commonly available.Therefore it is appropriate to schedule specific LiDAR acquisition flights.
In further analysis, if more detailed LiDAR data will be avail-able, it will be possible to produce a better 3D roof model and estimate an even more accurate solar potential of each building, and also high vegetation will be taken in account and its shadow introduced in solar radiation computation.

Figure 1 :
Figure 1: Web application Bologna Solar City

Figure 2 :•
Figure 2: LiDAR data coverage and Bologna district boundaries The laser scanner used by Blom CGR S.p.a. in the flight that took place in 2009 was a Optech ALTM 3033 with the following specification: Elevation points with orthometric height • District boundaries This dataset has a different CRS than LiDAR data: ED50 / UTM zone 32N (EPSG:23032) and orthometric height.In order to work in an accurate and coherent environment, all shapefiles were transformed into the same CRS as LiDAR data using ConvER3 GPS7, a free software from Regione Emilia-Romagna.It can use IGM (Italian Military Geographical Institute) grid to perform a rigorous transformation from different Coordinate Systems, including conversion from geoid to ellipsoid heights and viceversa.

Figure 3 :
Figure 3: Results from different Plane Surface Tolerance values

Figure 4 :
Figure 4: Different PST value with a curved surface

Figure 9 :
Figure 9: 3D Extracted roofs along with buildings footprints from pre-existing digital cartography

Figure 10 :
Figure 10: Building passed QA with little manual editing

Figure 11 :
Figure 11: 3D Viewer interface: details from the study area