This invention relates to a method for fitting three-dimensional data to a capped cylinder for use, for example, with a robotic tank inspection system, and more particularly to such a method that does not require iteration to find the best fit.
Liquefied petroleum gas (LPG) tanks are hazardous environments for inspection. LPG is cryogenically cooled into a liquid and stored in large (approximately 98 m diameter) tanks. These tanks must be periodically inspected, but pose significant hazards to the inspector. The temperatures are around −160 C and the fumes are explosive. The tank sits on a concrete slab so the floor is not accessible from the outside. Inspection must be done from inside.
Current methods to inspect the floor of LPG tanks are expensive and energy intensive [1]. The first step is to remove the LPG. The tank is then gradually warmed and a human inspector is 15 s sent in. The process of warming the tank to temperatures safe to send in an inspector can take up to two weeks due to the need to gradually warm the tank to avoid explosions. The inspection and cooling of the tank adds further time requirements. This leaves the tank out of commission for several weeks resulting in lost revenue. The heating and cooling is also expensive and uses significant energy.
A robotic inspection system will make the inspection system much more efficient. The tank will not need to be warmed as much as is necessary for a human inspector. This will shorten the down time, thus less revenue will be lost. This also means that there will be less energy necessary to cool the tank back down to operating temperatures. The inspection process itself can be done more quickly and thoroughly because the robots will move along prescribed paths more accurately than human inspectors. The risk is also diminished because robots are not sensitive to the same hazards as human inspectors.
Mobile robots often use range measuring sensors to determine the extents of their environment [2]. Advances in laser range finders have made them commonplace in mobile robot navigation systems. The data collected by these sensors need to be processed in order to create 3D maps of the environment.
Processing the data most often involves segmentation and then fitting each section to mathematical descriptions, mainly geometric primitives [3][6]. The mathematical descriptions of each geometric primitive make up the map of the surrounding environment which can then be used in movement or other interaction algorithms. The basic segmentation methods used by most current methods of segmentation can be broadly sorted into three categories: Clustering, Seed and Grow, and Hough Transform.
In clustering, properties of each data point are compared. Points with similar properties are grouped together through iteration. Usually these properties are similarity in position and orientation of the surface normal at that point [7]. Clusters are considered to be parts of separate surfaces when the properties in that cluster lie outside the user defined threshold for similarity.
Seed and Grow is similar to clustering in that it groups together points with similar properties, however to find these features, random “seeds” are selected and neighboring points are compared and collected until there is a large discrepancy. Once the neighboring points of the growth section fall outside of a threshold for similarity, another seed is chosen and the process repeats until all points are used. A final sweep combines homogeneous regions [6][7].
The Hough Transform relies on a voting procedure to find the most common geometric parameters for a given area. The true parameters of the geometric primitive composing the surface of interest should show up as the most common measurements, but may be distributed by noise. A histogram-based analysis of the parameters will show the dominant value, corresponding to the actual parameters of the geometric primitive. [8][9].
These methods are very flexible, but there is room for improvement when segmenting inspection environment data. Processing the data from range sensors into a map is often computationally expensive and time consuming. The resource demands are particularly large when the environment is completely unknown. Most methods using the above processes of segmentation assume nothing about the environment from which the data they are segmenting comes. This makes the algorithms very flexible because no external information is needed and they are not limiting what the fit outcome should be. The downside is that they all rely on iteration to find the best fit. Fortunately, the general structure of the inspection environment is known in the case of inspection robots. This information may be used to provide basic estimates of shape that reduce the need for computationally expensive iteration. This shape knowledge offers the potential to make the new algorithm significantly faster than the existing segmentation methods.
The method for mapping capped cylindrical environments according to an aspect of the invention includes collecting data from a distance range sensor within a capped cylinder. A group of N neighbors surrounding a point of interest at a radius R is collected and the normal at each sample data point is found from the N neighbors. The data is then segmented into groups and planes are fitted to the roof and floor and a circle is fitted to wall data. In a preferred embodiment, the sensor is a laser range sensor. A suitable method for finding the normal uses least squares plane filling. In this embodiment, the segmenting step including aligning coordinate frames at the cylinder and sensor locations. It is preferred that the fitting step uses a Gauss image.
In another preferred embodiment, the number of neighbors is at least three. It is also preferred that the method include the step of inputting constants based on fit requirements. Suitable constants are accuracy of desired fit; sensor resolution; encoder step resolution, estimate of cylinder height; and estimation of cylinder radius.
The algorithm implementing the method of the invention offers improved robustness, accuracy and speed over existing segmentation and fitting methods. These improved properties are achieved through the use of a Gauss map, three dimensional histogram and shape knowledge driven segmentation.
a and 2b are schematic illustrations of the mapping of an uncapped cylinder and a plane onto a Gaussian sphere.
a and 3b are simulated raw data and a histogram of Gaussian sphere normals.
Here we describe the algorithm designed to fit planar capped cylinders to data collected from a laser range finder. Consider a cylindrical tank 10 and laser range sensor 12 located on its wall as depicted in
The orientation of the sensor as it moves during data collection is measured by the amount of rotation in φ and θ. The distance from the sensor to the surface it is measuring, along the z axis, will be referred to as d. Some aspects of calculation are easier to accomplish in a 110 rectangular coordinate frame. In these instances the sensor rectangular coordinate frame designated by x,y,z can be derived from θ, φ and z. xyz and θφz are related by:
x=d sin φ cos θ
y=d sin φ sin θ
z=−d cos φ (1)
The best data collection process for this algorithm involves a nearest-neighbors method [10]. Sampled points are evenly distributed about a unit cube in the xyz coordinate frame. After the location of the sensor, θ and φ, and the distance to the surface, d, for each sample is collected, a group of N neighbors surrounding the point of interest at a radius R are also collected. This neighborhood is used to find the normal at each sampled data point using least squares plane fitting. Previous implementations of this method have determined N and R through trial and error [10]. A noise analysis was carried out to determine the best values and is discussed below.
There are three steps of segmentation which follow data collection. The first step of segmentation is locating the orientation of the Z axis in terms of θφz. Next the Z axis is aligned with the z axis through coordinate transforms. Once the two frames are aligned the data is segmented into groups corresponding to the floor, roof and walls. The fitting process consists of fitting planes to the roof and floor and a circle to the wall data. These three fits in addition to the rotation of the sensor origin through the transformation matrix defines the extents of the tank.
The process of aligning the coordinate frames then breaking the data into geometric primitives was chosen over fitting an arbitrarily rotated cylinder for several reasons. The first is to simplify the segmentation process. Traditional least squares fitting of a 3D rotated cylinder in space cannot handle the data corresponding to the capped ends. The excess data biases the result, ruining the cylinder fit. This suggests that the roof and floor data should be separated from the wall data. The segmentation methods listed above can accomplish this, but they are computationally expensive due to their iteration steps. Aligning the two coordinate frames places the cylinder (the XYZ frame) into a known arrangement. This allows us to remove the need for iteration, greatly reducing the computational load and speeding up the process. Another advantage of the rotation method is that non-uniformly distributed data does not cause a bias fit. The lack of bias is due to the location based method of segmenting the data points instead of the common property methods of segmentation. In addition, the walls are known to be vertical, thus we can constrain results to that orientation, reducing the degrees of freedom of the fitting algorithm and thus improving the accuracy of the fit.
While fitting data to cylinders has been done using multiple methods, all these processes have looked at data where the cylinder is modeled as infinite. The ends, or caps, either do not exist, as in the case of fitting to pipes in a building [5], or they have been processed separately [6]. One aspect that is unique to the method discussed herein is the integration of the caps in the cylinder fitting. This is achieved through the use of the Gauss image.
The Gauss image is formed through the placement of the ends of normals, collected from a surface, onto the origin. These normals then intersect a unit sphere, referred to as the Gaussian sphere, forming a Gauss map [4]. The resulting image when the entire surface is mapped is the Gaussian image. The Gaussian image of an uncapped cylinder is a circle on the Gaussian sphere as shown in
This algorithm is intended to work with cylinders composed of planar caps. The cylinder is then a combination of two planes and a cylinder. This means one can expect the Gaussian image to be a combination of the two images shown in
The coordinates of the points on the Gaussian sphere are collected and sorted into a 3D histogram in order to identify the Z axis. The histogram is composed of bins corresponding to the coordinates of the points on the Gaussian sphere. The histogram bin containing the most “votes” will correspond to the normal of the cap.
The cylinder is oriented in a known arrangement through the alignment of the sensor and global coordinate frames. The sensor z axis is aligned with
rotaxis=
and the angle about which to rotate is found from:
θrot=cos−1(
A rotation matrix using the rotaxis and θrot orients the cylinder along the z axis which is what allows for the non-iterative segmentation.
In order to easily fit the floor and roof planes and the cylinder to the walls, the data needs to be sectioned into the corresponding regions. Once the data has been rotated, the floor data can be identified as the data with z values below the clean floor height and with x and y values that are within the clean floor radius. The clean floor height and radius are calculated by modifying the 5th and 95th percentile values of the data according to the acceptable estimate error.
The isolation of the floor data is shown in
The same isolation process is completed for the roof, except the only values collected are those above a calculated clean roof height and within the clean radius. The wall data is collected by finding points that are less than the clean roof height, but above the clean floor height and are outside the clean radius.
Simple least squares fitting of each set of data results in the tank map. The floor and roof sets are fit to planes as shown in the top row of
The algorithm described above requires the designer to input five constants based on the fit requirements. These are 1) the accuracy of the desired fit, Acc, 2) the resolution of the distance sensor, σs, 3) the resolution of the encoder, σe, 4) a rough estimate of the height of the tank, Hrough, and 5) a rough estimate of the radius of the tank, Rrough. These inputs are used to calculate three variables: the neighborhood radius R, the number of neighbors N, and the total number of normals, Ntot which are used in the identification of
R is determined by the size of the tank and noise in the sensor. The neighborhood plane fit error will become relatively large if the R is too small and the normal fit will suffer. There will be a smoothing effect where the walls and corners will generate spurious normals if R is too large. The designer should use the largest R that will not result in excessive smoothing.
The number of neighbors, N, used in the normal calculation also affects the accuracy of the fit and the time to run. The absolute minimum number of neighbors is three in order to guarantee enough information for a plane fit.
More than three neighbors increases the accuracy as the noise is averaged out, but it takes more time to run. The optimal N is the minimum number of neighbors for which the algorithm will achieve the desired accuracy or three, whichever is larger. The desired accuracy, Acc, must be converted to the allowable angular error of the estimated tank floor normal and the actual tank floor normal, Err, to determine the N needed to achieve the desired accuracy as shown in equation 4,
An equation for N can be calculated from the sensor and encoder errors, σsand σe, the neighborhood radius, R, the estimated height and radius, Hrough and Rrough, the distance of the sensor from the wall as determined by the geometry of the sensor platform, s, and the acceptable angular error Err using equation 5.
R and N can be derived from the user inputs as described above, or the user can modify them directly. These values are derived using equations that optimize the fit at the expense of time. If time is more critical to the application, the designer can tweak the values of R and N down to decrease the time. If the sensor is ideally located, these numbers could be very small. If the sensor is located in a difficult position, such as near a corner, there may not be much the designer can do without losing accuracy.
Although not an input, the number of normals to find before moving on to the segmentation step, Ntot, is important to the accuracy of the fit. This forms another tradeoff between time and accuracy. Measurements of about 50 normals have been found to work in most cases, but significantly higher values (upwards of 200) are needed when the sensor is placed in the extreme corners of the tank. As a default, 218 normals are used in order to guarantee a good fit. The designer can lower this value to achieve better processing time if the sensor is not badly located.
The algorithm disclosed herein was developed and tested using simulated data. The simulation was able to produce data with varying origin locations, tank heights and radiuses, cylindrical rotation and noise values. The input characteristics gave a ground plane against which the fit may be tested. The algorithm is able to fit a mathematical model within the specified accuracy given variation in all parameters.
The contents of all of the references cited herein are incorporated into this application by reference.