CITY STRUCTURE ANALYSIS ON QUICKBIRD IMAGERY BY MULTISCALE RADON TRANSFORMATION

The radon transformation is a proven image analysis technique in medical image processing, although there are some applications with satellite images. Among these applications maybe the most interesting one is based on the feature that this mathematical transformation is capable of detecting roads and road junctions. This option makes the transformation similar to the well-known Hough transformation, but radon has better performance parameters. Roads are ground objects having linear element borderlines, which can be extracted from radon analysis. Since roads and streets are crucial objects in urban environment analysis, they are excellent sources to get information about the city structure. The current paper introduces the short mathematical background behind the technique, and then a multiscale approach is described. A Quickbird image captured about Budapest was used to create a test data set with different scales, where the radon transform is applied. The different scale image layers show different details about the road structure, therefore about the city architectural structure. The result of the analysis can be presented in GIS database, thus the extracted information is able to be reused in further sophisticated applications. 1. INTRODUCATION The modern high resolution satellite sensors are able to generate data about the Earth’s surface in an extreme amount. The data analyzing techniques must therefore cope with the data amount and have to extract relevant information. The Quickbird, Ikonos, WorldView and other satellite imaging tools monitor the continents with high geometrical resolution and repetition rate. In the current paper the urban pattern description is aimed to be supported by image processing techniques, taken from other imaging field. The author performs a trial to elaborate a technology to extract land cover information focusing on urban areas. The algorithm is “under construction”, but the actual results seem to be interesting for the photogrammeters, too. 2. THE RADON TRANSFORMATION The radon transformation is a mathematical projection of an image along several specified directions. The transformation computes a line integral from multiple sources along parallel paths (the so-called beams) of an f(x,y) two-dimensional image matrix. Fig. 1 shows how an arbitrary object can be projected by parallel beams. The positioning of the parallel paths of the beams is given by an angle (theta), similarly like a rotation is defined. All integrations happen along the paths. If the image and its content is represented by the function f(x,y), where x and y are the image coordinates, the integration is executed along a transformed axis. The transformation is a rotation, which can be easily described by a rotation equation: ( ) ∫ ( ) (1) Figure 1. Schema of the projection with parallel paths (beam) In the above equation the two-variable function requires the rotation, as given


INTRODUCATION
The modern high resolution satellite sensors are able to generate data about the Earth's surface in an extreme amount.The data analyzing techniques must therefore cope with the data amount and have to extract relevant information.The Quickbird, Ikonos, WorldView and other satellite imaging tools monitor the continents with high geometrical resolution and repetition rate.In the current paper the urban pattern description is aimed to be supported by image processing techniques, taken from other imaging field.The author performs a trial to elaborate a technology to extract land cover information focusing on urban areas.The algorithm is "under construction", but the actual results seem to be interesting for the photogrammeters, too.

THE RADON TRANSFORMATION
The radon transformation is a mathematical projection of an image along several specified directions.The transformation computes a line integral from multiple sources along parallel paths (the so-called beams) of an f(x,y) two-dimensional image matrix.Fig. 1 shows how an arbitrary object can be projected by parallel beams.The positioning of the parallel paths of the beams is given by an angle (theta), similarly like a rotation is defined.All integrations happen along the paths.If the image and its content is represented by the function f(x,y), where x and y are the image coordinates, the integration is executed along a transformed axis.The transformation is a rotation, which can be easily described by a rotation equation: Figure 1.Schema of the projection with parallel paths (beam) In the above equation the two-variable function requires the rotation, as given The inverse of this rotation is also necessary, but luckily it's very easy, because of the orthogonality of the transformation: As these arguments are written into the function f(x,y), the integration formula can be present:  Remarkable applications of the inverse radon transformation can be found in object reconstruction, in computer vision, and in imaging of computer tomography.One can extract also a subset of the accumulator image; in this way a selected set of lines are reconstructed.Following this strategy, the most significant (the strongest) lines can be drawn into an image, which is useful in city structure analysis, too.

TEST DATA
The applied test data is acquired by the Quickbird satellite.The sensor acquires images about Earth's surface in panchromatic and multispectral range.The panchromatic sensitivity has its maximum near to the wavelength of 725 nm.The geometric resolution is about 61 cm in nadir direction.The radiometric resolution has been defined by an 11-bit analogue-digital conversion.

IMAGE ANALYSIS
Prior to the radon transformation several steps were necessary.The first was the resampling, where 1m, 3m, 5m and 10 m resolution images were derived.The resampling used the bicubic method in order to get smooth intensity values but less contrast loss.The input images for radon transformation were computed by a very simple way: an empirically defined threshold was applied, then the achieved binary image was filtered by a circular average filter with a size of 3 to 5 pixels.This procedure roughly extracts the building roof pixels.Of course more sophisticated edge detection, like Canny or Deriche-filters can produced more precise input for the transformation.The radon transformation has an angular resolution of 1 degree, the range is from 0 to 179 degree.This interval covers the possible lines without angular redundancy; skipping the 180-360 degree interval accelerates the work.A very convenient feature of the radon transformation is its rotation invariance.If a structured image would be rotated and radon transformed, the achieved result will differ only in x-shift.Therefore the line crosses will produce the same result image, independently of their orientation in the original image.Furthermore, the radon image shows correlation with the details in the input images: the higher the geometric resolution is, the more scattering can be observed.The reason is that the geometric resolution controls the pixel size, therefore by continuously reducing the resolution, the pixel size increases and the image details are getting lost.After losing the significant fine details, the ability of detecting lines, therefore the streets and hence the city structure will be lost.The high resolution satellite images contain sufficient geometric details to detect the streets.The line extraction, i.e. the reconstruction must be done by selecting the highest value combinations.To achieve a small smoothing in the radon image, a median filter was applied.As an example, the most relevant 30 lines are reconstructed and projected onto the Maglód-area image.

CONCLUSION
The current research focuses on the geometric image transformation, called radon transformation.The algorithm converts the binary image content into a special polar coordinate system, where the possible line elements are represented by the orientation angle and their distance to the image center.If one can extract the significant image element, e.g. by thresholding a gray scale image into a binary one, the radon transformed result contains the accumulated lines.By selecting the highest peaks in the accumulator image, the most relevant directions can be detected.The experiment has proven that city subsets with different densities have different radon counterpart, where the city density strongly correlates with the radon image.The test also proved that the method depends on the information content of the image, i.e. on the geometric resolution or on the pixel size.The applied Quickbird imagery has the adequate information content to extract such city structure information.The future application of the derived achievements could be a segmentation technique of urban areas, where suitable sized running window cuts off the test fields and the connected radon transform gives back the most relevant building directions.An adequate categorization will show the architectural structure of the working area.

Figure 2 .
Figure 2. The rotated coordinate system in radon transform ( ) ∫ ( ) (4) This theoretical formula is practically realized by summation with a limited resolution along y'.The angular resolution (along theta) is discrete, of course.The radon transformation converts the intensity image into a kind of accumulator image where the axes represent the x' and theta coordinates.Examples will be shown in Chapter 4.

Figure 3 .
Figure 3.An arbitrary angle slice about a rectangle object The inverse radon transformation computes the intensity values from the accumulated values considering the rotation formula.Remarkable applications of the inverse radon transformation can be found in object reconstruction, in computer vision, and in imaging of computer tomography.One can extract also a subset of the accumulator image; in this way a selected set of lines are reconstructed.Following this strategy, the most significant (the strongest) lines can be drawn into an image, which is useful in city structure analysis, too.

Figure 6 .
Figure 6.Radon images of the TA1 sample with different geometric resolution Focusing on only one test area, different resolution sampling can executed.By this multiscale or multi resolution analysis the sensitivity of the method can be studied.Fig. 6 presents the same area with 1-3-5-10 m resolution, respectively.The scales on the left of the diagrams have different ranges.This is because of the modified number of lines in the input binary image.The vertical resolution of the radon image equals to the number of rows of the input image.The horizontal size of all radon images are the same: they are set to the angular resolution.Furthermore, the radon image shows correlation with the details in the input images: the higher the geometric resolution is, the more scattering can be observed.The reason is that the geometric resolution controls the pixel size, therefore by continuously reducing the resolution, the pixel size increases and the image details are getting lost.After losing the significant fine details, the ability of detecting lines, therefore the streets and hence the city structure will be lost.The high resolution satellite images contain sufficient geometric details to detect the streets.The line extraction, i.e. the reconstruction must be done by selecting the highest value combinations.To achieve a small smoothing in the radon image, a median filter was applied.As an example, the most relevant 30 lines are reconstructed and projected onto the Maglód-area image.

Figure 7 .
Figure 7. Extracted lines projected onto the 3 m resolution image The image processing was implemented in Matlab 2010b with Image Processing Toolbox on an Intel Mobile Core2 Duo 2.53GHz computer.The software can use the acceleration of the graphics card (NVidia CUDA), but no extra GPU optimization was required.The applied processing functions needs several seconds computing time.

Table 1 .
Description of the satellite samples The imagery was captured inMay 2002, August 2004 and  October 2005.Only the panchromatic images were used in the current research.The city structure analysis has been started by selecting five land cover categories in the agglomeration of Budapest.These categories are described in Table1.