ESTIMATION OF TREE STEM ATTRIBUTES USING TERRESTRIAL PHOTOGRAMMETRY

The objective of this work was to create a method to measure stem attributes of standing trees on field plots in the forest using terrestrial photogrammetry. The primary attributes of interest are the position and the diameter at breast height (DBH). The developed method creates point clouds from images from three or more calibrated cameras attached to a calibrated rig. SIFT features in multiple images in combination with epipolar line filtering are used to make high quality matching in images with large amounts of similar features and many occlusions. After projection of the point cloud to a simulated ground plane, RANSAC filtering is applied, followed by circle fitting to the remaining points. To evaluate the algorithm, a camera rig of five Canon digital system cameras with a baseline of 123 cm and up to 40 cm offset in height was constructed. The rig was used in a field campaign at the Remningstorp forest test area in southern Sweden. Ground truth was collected manually by surveying and manual measurements. Initial results show estimated tree stem diameters within 10% of the ground truth. This suggest that terrestrial photogrammetry is a viable method to measure tree stem diameters on circular field plots.


Background and Aim
The forestry industry requires information, such as species, stem diameter, volume and height, about trees in the forests for management planning.Traditionally, the information has been collected by manual inventory using calipers to measure stem diameters.The inventory in Sweden is usually done on circular field plots with a radius of 10 meters.From this kind of samples larger stocks of forests are estimated.The process of manually measuring the trees is labor intensive, and the planners want as much and as reliable information as possible.
In this paper, a method to estimate tree stem diameter and position using terrestrial photogrammetry and image analysis is proposed.

Related work
There are only a few studies to be found on estimation of stem attributes using terrestrial stereophotogrammetry.The work by Reidelstürz (1997) was performed on analog images and concludes that manual measurements in analog images would be unfeasible due to long processing times.More recently, Hapca et al. (2007Hapca et al. ( , 2008) ) present methods for stem profiling and measurement based on digital photogrammetry.However, this method places several constraints on the imaging process.The cameras must be accurately leveled and the images must be taken from two angles, separated by exactly 90 .Furthermore, the scene must be modified to include a marking of a known distance.The paper by Abd-Elrahman et al. (2010) presents a method for inventory of trees in urban areas using images taken by the general public.The presented method is based on manual measurements in digital images and is evaluated on free-standing trees in a number of parks.The analysis time for a set of 6-12 images is cited as ⇤ Corresponding author 4-8 hours.The report by Larsson et al. (2011) studies the related problem of sensor technology to be mounted on harvesters to collect forestry data during harvesting.The sensors include a pair of video cameras, a 2D line scanner (SICK LMS-511), and a short range flash sensor (PMDTec Camcube).

MATERIALS AND METHOD
2.1 Camera rig 2.1.1Hardware A camera rig consisting of five digital system cameras with a total baseline of 123 cm and up to 40 cm offset in height was built (Figure 1).Three Canon EOS 7D equipped with Sigma EX20/1.8DG wide angle lenses and two Canon EOS 40D equipped with Canon EF-s 17-55/2.8zoom lenses, locked at 17 mm, were mounted to the rig.The cameras were synchronized by a cable trigger to ensure the images were taken at the same time.A magnetic field compass was attached to the rig to make it easier to evenly distribute the exposures in a circle.
Five cameras were used for two reasons.Firstly, offset-mounted additional cameras reduce the risk of false matches due to ambiguities in epipolar matching between only two cameras.Secondly, the larger number of cameras reduces the risk of occlusion problems, in an environment where branches and leaves at short range from the cameras are common.
The total weight of the rig was about 13.5 kg.The rig can be used with a common camera tripod.However, this would be cumbersome to carry around and adjust in the terrain.Instead, Easyrig 2.51 , a portable camera support system originally developed for TV cameras, was preferred (Figure 1).

Calibration
To investigate the stability of the camera and lens before they were mounted on the rig, each of the cameras was calibrated multiple times using the PhotoModeler camera calibration sheet 2 and in-house camera calibration software (Börlin et al., 2004).For additional stability, the calibration sheet (Figure 2) was printed on an aluminium sheet with attached stringers.Over 20 calibration images were used for each camera calibration to cover different camera rotations, camera-to-object distances, and all of the CCD.The cameras were calibrated four to ten times each during a period of two months.The variation in focal length estimation during this period was below 0.03 mm.
The camera rig was calibrated by bundle adjustment of one exposure of the same calibration sheet.The calibration was repeated on a regular basis during the campaign to reduce any handling effects.

Field work
2.2.1 Data acquisition protocol Within each circular field plot, images were acquired from five positions.At every position, exposures were taken in each of 12 directions separated by approximately 30 degrees as shown in Figure 3. Two sets of images were taken in each direction to increase redundancy and to guard against problems with the exposure.The camera rig was calibrated before and after every plot using the procedure described 2 www.photomodeler.com in Section 2.1.2.Each exposure was furthermore recorded in a written protocol.

Campaign
In August 2011, a field campaign was undertaken to the Remningstorp test area outside Skövde in southern Sweden (N58.46 , E013.65 ).The forest at Remningstorp is semi-boreal with mostly Scotch Pine, Norway Spruce and some broad-leaved trees.
During three days of field work, 25 field plots with 20 m radius were photographed using the five-camera rig.An example of a single five camera exposure is shown in Figure 4.The weather during the campaign was mostly sunny with scattered clouds, with some period of broken clouds and light rain.In total, about 25 000 images summing up to 175 GB of data was collected.
At each field plot, a stick in the ground marked the center of the plot.Images were acquired as close as possible to the center.Occasionally, obstacles like big stones, a stream, or a tree close to the center caused the camera station to be moved to a usable position.The four satellite camera stations were found using the compass and stepped out by someone in the team.Occasional reference measurements by a tape measure were made at some places, indicating that the satellite distance was usually within one meter of the correct value.Approximately 20-30 minutes were spent on each plot.This time covers image acquisition including rig pre-and post-calibration, transportation to and from the closest road as well as loading and unloading of the equipment from the car.

Ground truth
For validation, the stem diameters at breast height (DBH) of all trees in the plots were manually measured separately using a caliper.Furthermore, tree positions within each plot were determined by a land surveyor using a total station.Additionally, the plots were laser scanned from the plot center with a Leica ScanStation C10.

Algorithms
From each five-tuple of images acquired by the rig, a point cloud is constructed from feature points detected by SIFT (Lowe, 2004).Epipolar constraints (Hartley and Zisserman, 2003) are used to minimize false matches.Figure 4 shows a set of images from the camera rig.In the generated point clouds, some vertical clusters representing tree stems are obvious (Figure 5).
The points are projected on a ground plane, simulated by the x z plane in the rig coordinate system.A rough estimation of the stem diameter at breast height is acquired by circle fitting using RANSAC (Fischler and Bolles, 1981), the algebraic fit method by Gander et al. (1994) and a geometric fit method using Gauss-Newton (Gander et al., 1994;Nocedal and Wright, 2006).
3.1.1Point cloud generation algorithm using N cameras from a single viewpoint Given detected SIFT features in each image, a point cloud is generated by the following algorithm: 1.For each feature point i in the reference camera r (typically the center camera) 1.1 Create a set of selected points S, initially containing the single point i from camera r.
1.2 For each other camera j 6 = r 1.2.1 Create a set of candidate matches Cj, initially containing all feature points in image j.
1.3 Select image r as the active image.
1.4 Repeat while there is at least one non-empty set Cj of candidate matches 1.4.1 For each non-empty set Cj 1.4.1.1Filter the candidate matches Cj by the selected point p in the active image.A candidate point q 2 Cj is removed if • the descriptors of q and p are too dissimilar, or • the distance from q to the epipolar line of p is too large, or • the q and p frame sizes differ too much, after compensation for the expected difference due to potentially differing image resolutions.1.4.2If there exists a set Cj consisting of a single point k, the matching in image j is unique: 1.4.2.1 Add the point k in image j to the set S of selected points.1.4.2.2 Clear the set Cj as the single remaining candidate has been selected.1.5 If points were selected in at least three images, calculate a 3D point by forward intersection (McGlone et al., 2004).
An example of a generated 3D point cloud is shown in Figure 5.

Tree stem diameter measurements
Given the 3D point cloud, tree diameters are calculated by the following algorithm: 1. Project the 3D point cloud to the x z ground plane (Figure 6).
2. Segment the point cloud to find groups of points corresponding to a single tree.This is currently done manually.
3. Fit a circle to each segmented point cloud using RANSAC.
4. Fine-tune the circle to the remaining point using a Gauss-Newton algorithm.

Data selection
The images acquired from the central camera station of one field plot have been analyzed.The chosen plot is typical in the sense of having a variation in tree size and includes parts with only low undergrowth and other parts with brushwood.Spruce, pine, birch and other leaves are represented.As an initial data set, a total of twelve trees were chosen as the most clearly visible trees in the 3D point clouds.

PRELIMINARY RESULTS
An initial implementation of the algorithm of Section 3.1 has been applied to images from the center camera station from one of the field plots.Figure 4 shows an image set with matched feature points marked.The green lines are epipolar lines for one of the points.The generated point clouds are sparse due to high values of the filtering thresholds in Step 1.4.1.1 of the algorithm, but enough points exist to make some trees easy to detect, see Figure 5.
In Figure 6, the points in Figure 5   has been manually indicated.In Figure 7, the indicated points have been filtered by RANSAC and circles have been fitted to the filtered points using three methods: The direct algebraic method by Gander et al. (1994) and two geometric fit methods using the Gauss-Newton algorithm (Gander et al., 1994;Nocedal and Wright, 2006).The unweighted geometric fit methods assumed isotropic point errors whereas the weighted geometric fit method assumed anisotropic point errors calculated by forward intersection.
A set of stem radii, estimated by the above algorithms, are presented in Table 1.All stems were positioned at a maximum distance of approximately seven meters from the rig.Ground truth from caliper measurement is also displayed.The method by Gander performed the best, with most of the results within 10% of the ground truth.Unweighted Gauss-Newton was almost as good, but the weighted Gauss-Newton overestimates most of the radii.Observe that tree number 11078 is underestimated in image set 19 and overestimated in image set 15. Manual inspection of matched points in the images indicates multiple points on the ground.

CONCLUSIONS AND FUTURE WORK
These preliminary results indicate that this method is viable for estimation of stem attributes in circular field plots.The method will be further developed and validated.
The current method appears to be sensitive to outliers remaining after the RANSAC filtering.This could be handled by levelling of the rig to obtain a better approximation of the ground plane followed by automatic exclusion of points too close to the ground.Circle fitting at different levels and/or fitting of a conical cylinder to the 3D point cloud could be used to handle the tapering of the stem.Furthermore, the RANSAC algorithm will be extended to handle anisotropic point errors.Additionally, the rig calibration will be improved to include multiple calibration images acquired at different distances to the calibration object.
One shortcoming in the current code is that small stems are hard to measure due to too few detected features.A possible improvement would be to co-register the point clouds from multiple exposures and camera positions in a field plot.The co-registration of the point clouds should also make it possible to create a tree map over the test site, which could be useful for different applications.Furthermore, comparisons of the stem diameter estimates will be performed with point clouds produced by a terrestrial laser scanner placed at the plot center.

Figure 1 :
Figure 1: The five-camera rig in the field, carried with Easyrig 2.5.

Figure 2 :
Figure 2: The PhotoModeler camera calibration sheet printed on a sheet of aluminium.Stringers were attached to the back of the sheet for additional stability.

Figure 3 :
Figure 3: The image acquisition protocol for each field plot consisted of a central camera station and four satellite stations.The satellite stations were approximately 10 m from the center in the four cardinal directions.At each station, images were taken in 12 directions separated by approximately 30 degrees.

Figure 5 :
Figure 5: Part of a 3D point cloud.The trunk of a tree is clearly seen.The axis unit is meters.1.4.2.3 Select image j as the active image.1.4.3 otherwise if there are non-empty sets left, the matching is ambiguous in the remaining image(s): 1.4.3.1 Abandon the remaining images by clearing the remaining non-empty sets Cj.
Figure 4: A set of images from the five-camera rig.The relative positioning of the images indicate the relative position of the camera within the rig.The green lines are the epipolar lines for one point.

Figure 6 :
Figure 6: The point cloud of Figure 5 projected on the x z plane.The points in the red rectangle correspond to points on the tree trunk.The axis unit is meters.

Figure 7 :
Figure7: Three circles fitted by different methods to the points selected in Figure6after RANSAC filtering.Each point is indicated by its 1-sigma uncertainty ellipse.Some remaining outliers are seen, e.g.inside the circles near (-1.05,4.9)or outside the circles near (-0.75,4.85).The blue circle is fitted using the direct algebraic method.The green circle was fitted by the Gauss-Newton algorithm assuming isotropic error.The red circle was fitted by the same algorithm assuming the errors indicated by the blue ellipses.The axis unit is meters.