Detecting Dead Standing Eucalypt Trees from Voxelised Full-Waveform Lidar Using Multi-Scale 3D-Windows for Tackling Height and Size Variations

In southern Australia, many native mammals and birds rely on hollows for sheltering, while hollows are more likely to exist on dead trees. Therefore, detection of dead trees could be useful in managing biodiversity. Detecting dead standing (snags) versus dead fallen trees (Coarse Woody Debris—CWD) is a very different task from a classification perspective. This study focuses on improving detection of dead standing eucalypt trees from full-waveform LiDAR. Eucalypt trees have irregular shapes making delineation of them challenging. Additionally, since the study area is a native forest, trees significantly vary in terms of height, density and size. Therefore, we need methods that will be resistant to those challenges. Previous study showed that detection of dead standing trees without tree delineation is possible. This was achieved by using single size 3D-windows to extract structural features from voxelised full-waveform LiDAR and characterise dead (positive samples) and live (negative samples) trees for training a classifier. This paper adds on by proposing the usage of multi-scale 3D-windows for tackling height and size variations of trees. Both the single 3D-windows approach and the new multi-scale 3D-windows approach were implemented for comparison purposes. The accuracy of the results was calculated using the precision and recall parameters and it was proven that the multi-scale 3D-windows approach performs better than the single size 3D-windows approach. This open ups possibilities for applying the proposed approach on other native forest related applications.


Introduction
Tree deaths are significant in a circular ecology, since they provide resources to many organisms [1] and habitat to many mammals and birds [2].Zielewska-Büttner et al. [3] showed that the abundance of dead standing trees (snags) was the most important predictor of woodpecker habitat sections, stressing this way the importance of dead wood.At the same time, decaying wood is decreased at managed forests and this is a thread for organisms whose lives depends on dead wood [4].According to Rose et al. [5] it is less expensive to preserve natural habitant than to regenerate it.There is, therefore, a need for understanding and maintaining the natural distribution of dead wood across a forest.
In Southern Australia, river red gums (Eucalyptus camaldulensis) grow along the banks of Murray river and its floodplains.Floods are key ecological components, while extinctions of native species could be driven by anthropogenic factors [6].After the construction of Lake Hume in 1934, the flow of the Murray river is regulated [7].As a result the annual flow of the river have been decreased [8] and fewer, shorter floods occur annually [9].George [10] identified increased stress and health decline of red gums, as well as reduction in seeds production, which threatens decline of the population.
Furthermore, tree hollows play a substantial role in preserving biodiversity in Southern Australia [2,11], because most arboreal mammals and numerous native protected bird species rely on them for sheltering [12].Nevertheless, in Australia there are no hollow creators like the woodpeckers that exists in the northern hemisphere.For that reason, it takes hundreds of years for a hollow to be formed [13] by insect and fungal attacks when access points are provided through damage caused by wind, storms.Furthermore, Gibbons et al. [14] claims that dead trees or trees in poor physiological condition are more likely to contain hollows.In parallel, Lindenmayer and Wood [15] and Goldingay [16] predicted that in the near future shortage of hollows for colonisation will exist due to anthropogenic factors.For these reasons, this study aims to contribute in managing dead wood in Southern Australia by improving automated detection of snags.
This study focuses on detecting dead standing eucalypt trees from a Southern Australian forest that has been influenced by the reduced flow and floods of Murray river.Even though detection of both fallen [17] and standing [18] dead trees is important for managing dead wood in Southern Australia and preserving biodiversity, they are addressed differently from a classification perspective.Fallen trees are detected by identifying line-like features on the Digital Terrain Model (DTM) that is created from LiDAR point cloud [19,20].In respect to dead standing trees, the following features could contribute into detecting them: their light reflectance since they absorb more green light [21] and their shape since they are less leafy and more likely to have broken branches [22].Additionally, Yao et al. [22] and Shendryk et al. [23] performed tree delineation before classifying standing trees as dead or alive.
LiDAR is extremely useful in forestry because the laser can penetrate the forest canopy through the gaps between branches and leaves.Therefore, significant structural information about forest structure at tree level are collected.LiDAR systems used to record only discrete peak returns (discrete LiDAR) that usually corresponds to large branches and the ground.In these discrete LiDAR systems, a minimum distance between two recorded returns existed; for the Leica ALS50 sensor there has to be at least a 2.7 m gap between two recorded returns [24].Over the years and with the technological advances, LiDAR systems become able to digitise and record the entire backscattered signal.The backscattered signal is digitised and stored into a number of waveform samples equally spaced (e.g., 15-30 cm vertical resolution).The intensity of each waveform sample corresponds to equal pulse width since the signal is digitised at equal space time intervals and, therefore, the time distance between every two coherent waveform samples is constant.This is explained in the file format specifications of LAS1.3 [25], LAS1.4 [26] and Pulsewaves [27].In 2006-2009, finding peak points from the waveform data [28,29] Reitberger et al. [30] Chauve et al. [31] attracted the interest of the scientific community.Scientists were able to find additional returns that were not acquired by the discrete systems due to the minimum gap that had to exist between two returns [28].Reducing the waveform data to discrete returns is easier to handle the increased amount of information recorder and work with existing work-flows.Nowadays, many sensors acquire only waveform data (e.g., Trimble AX60) and discrete data are produced by analysing the waveforms and extracting peak points at post-processing.Therefore, the terms "extraction of peak returns from full-waveform (FW) LiDAR" and "discrete LiDAR" can sometimes be used interchangeably.Even though, Anderson et al. [32] proved that FW LiDAR worth the extra processing and can estimate forest related parameters better, there are still questions to be answered; with the increased acquired pulse density and the advancement of technologies does the waveform still worth it?In addition, even though extraction of peak points from FW LiDAR could be identical to the delivered discrete data (tested by the authors), there are alternative ways of interpreting the waveform data, e.g., classification of waveforms according to their shape [33] and voxelisation [34].
Voxelisation is the process of inserting either the waveform samples or the discrete points into a 3D regular grid.Afterwards this voxelised data are used to derive terrain, canopy and other tree related metrics.The concepts of voxelisation (as explained below in Section 2.2.1 or with similar interpretations) have been used in forestry for handling both discrete [35][36][37][38][39][40] and waveform [41][42][43] data.In comparison to the discrete data, the waveform data are more likely to contain noise.Nevertheless, pulse density depends on the speed of the flight and the scanning pattern, while estimation of forest parameters is pulse density dependant [44].The intensity of each waveform sample corresponds to equal pulse width [25][26][27].They are, therefore, comparable to each other and the intensity values of each voxels can be normalised, overcoming the uneven density of LiDAR footprints [45].This study focuses on the algorithm that tackle height variations while working with eucalypt native forests and, therefore, comparison of the results between discrete and FW LiDAR data is not conducted.It uses voxelisation for interpreting the FW LiDAR data, but it could have voxelised the discrete LiDAR instead and applied the same methodology.
Classifications at tree level, while working with a native forest is a challenge since tree delineation is usually performed before health assessment [22,23].Tree delineation can be achieved by firstly detecting local maxima from the Canopy Height Model (CHM) and then segmenting CHM into individual trees with the watershed algorithm [46].The introduction of the marker controlled watershed algorithm improved the results [47] and further improvements were made by including structural information of the tree trunks and the under-storey layers of the canopy [48].Bottom-to-top delineation was proposed by Lu et al. [49] for segmenting deciduous trees from data collected during the leaf-off season.Similarly, Shendryk et al. [50] published an interesting bottom-to-top red gum delineation algorithm [22].Once trees are delineated, they can be classified as either dead or alive.
Previous work used multi-scale 2D analysis on Digital Elevation Models (DEM) or Canopy Height Models (CHMs) for detecting tree tops and delineating trees [51].Jing et al. [47] also used a multi-scale 2D approach on the CHM of dense points cloud (45 points/m 2 ) for tree delineation: at first Jing et al. [47] used scale analysis for determining dominant tree size.Then, they produced segmentation maps at multiple scales using the marker-controlled watershed algorithm and, finally, fused the multi-layered segments.Hu et al. [52] improved this approach by using Gaussian analysis to determine whether a segment consists of multiple trees.
This paper attempts to address the limitations of working with eucalypt trees in a native Australian forest.Detection of dead standing eucalypt trees from full-waveform LiDAR without tree delineation has been proposed before and it was shown that it performs better than a random prediction [53].In image processing and computer vision, scientists try to identify if objects, like faces, exist within a 2D image and detect them by extracting features [54].Similarly, Miltiadou et al. [53] extract 3D structural features from local areas around dead trees (positive samples) and live trees (negative samples) using single size 3D-windows and trains an object detection classifier.For tackling height and size variations of trees, this study adds on existing knowledge by introducing the usage of multi-scale 3D-windows.The classifier creates three probabilistic fields using three different sizes of 3D-windows and then merge the results, before proceeding to thresholding, filtering and assignment of predicted locations of dead standing trees.Using cross validation and comparison with the single size 3D-windows approach, it was proven that in comparison to a single-size 3D-windows approach, the multi-scale 3D-windows methodology improves prediction.

Study Area
Detection of dead standing eucalypt trees can be used as a proxy for managing biodiversity, since are more likely to have hollows for sheltering animals [12,15].It worth clarifying though that having a certain number of standing dead trees is important for biodiversity, but a whole forest of dead trees wouldn't make for a very biodiverse landscape.The study area is in south-eastern Australia, as shown in Figure 1 and has mild elevation differences (i.e., the area is not mountainous).According to "Forestry Corporation of NSW", Australia, it consists of 91.59% river red gum (Eucalyptus camaldulensis).The rest of the trees are black box (Eucalyptus largiflorens) and wattle group (Acacia spp.).Trees in this native forest are highly variable in terms of distribution, density and health state, because their regeneration depends extremely on floods [55].Extraction of forest related parameters from a native eucalypt forest is a challenge for multiple reasons.Eucalyptus camaldulensis tree height reaches up to around 40 m and they have multiple trunk splits.According to Wilson [56], their structural complexity makes interpretation of the acquired remotely sensed data difficult.To be more precise: (1) The sizes, heights and shapes of the trees significantly vary.It is, therefore, difficult to fit the features extracted about the dead trees into a single statistical model.(2) Eucalypt trees have multiple trunk splits and, therefore, each tree appears to have multiple local maxima on the CHM-tested by the authors.(3) The density of the red gum forest along the Murray river significantly varies since their regeneration depends on floods [55].For that reason a single type of tree delineation algorithm is not applicable for the entire study area.Figure 2 demonstrates the complexity of detecting dead standing eucalypt trees from the acquired Airborne Laser Scanning data.It is shown that individual trees are not distinguishable from the Canopy Height Model (CHM).There is also no consistency between the shapes of dead trees.There is, therefore, a need for advancing algorithms exploiting structural information for tackling height variations and variance forest density in native forests.

Acquired Full-Waveform LiDAR
Airborne LiDAR data were acquired by a Trimble AX60 LiDAR instrument.The data were acquired from March 6th to March 31st of 2015.The entire area of interest was scanned (Figure 1), resulting in around 4 TB of LiDAR data.The acquired data were subdivided into 206 flightlines and stored into LAS1.3file format [25].There was also a 30% swath overlap.The requested point spacing along and across the track was 0.48 m and the average spacing was 4.3 footprints per m 2 , including the denser areas due to flight overlapping.It is worth mentioning that due the size of the eucalypt trees (tall trees reach up to around 40 m) and the requested spacing, the system failed to acquire information about the trunks of the trees (see Figure 3).As a result, the state-of-the-art bottom-to-up tree delineation approaches that start by detecting the tree trunk and require denser footprints of LiDAR pulses could not be applied.LiDAR systems emit laser beams and collect information by measuring the round trip time of the pulse.The waveform digitiser records the intensity returned at equal time space intervals.Within the LAS1.3 [25] file, each waveform is attached to an anchor point that is calculated using three parameters: (1) a discrete return record (X record , Y record , Z record ), (2) a scale factor (X scale , Y scale , Z scale ) and (3) an offset (X o f f set , Y o f f set ) that defines the distance between the anchor point and the associated discrete point.According to [24] the recording of the waveforms starts a bit earlier than the first peak return recorded.This is achieved using a buffer since the recording is triggered at the first peak return.The anchor point is calculated as follows: Each waveform sample is equal to the laser intensity recorded during a time interval that is constant during digitisation [24].For each waveform recorded, there is an array containing the waveform sample intensities recorded.An example of the estimated average location of waveform samples is depicted in Figure 4.The average locations are calculated using the vector (X t , Y t , Zt), which indicates direction and distance between coherent waveform samples:

Field Data
The field data, provided by "Interpine Group Ltd.", New Zealand, were collected in July 2015 during the winter season of Australia.Field data measurements were taken from 33 plots with radius 35.68 m and area 0.4 ha.The field plots were randomly scattered inside the scanned area, which is coloured with red in Figure 1.On these randomly distributed plots, individual measurements of 2386 trees (of which 260 are dead standing trees) are stored.For each tree the bearing from the north and the distance from the centre of the plot are given.By knowing the position of the centre of the plot, the northing and easting co-ordinates of each tree are calculated.The locations of the field plots are depicted in Figure 1.Furthermore, they are spatially well-distributed and represent the entire study area; the minimum distance between two field plots was 3.1 km and the average distance between a plot and its nearest one was 7.8 km.
Due to the flights overlapping at the edges of the swath, there are 17 plots that were scanned twice from different angles and opposite sides.This resulted in significant differences between the overlapping acquired data and are, therefore, considered as different testing/training field data.This resulted in 50 unique associations between field plots and flightlines.These 50 associations are referred as field plots in the evaluation and methodology sections.In this study the northing and easting location of each dead and live tree was used to train the classifier.The locations of the trees within a sample field plot are given in Figure 2; red stars represent locations of live trees and the yellow stars represent locations of dead trees.
In this study, it is important to understand the variance of tree heights.In the field, only a sample of tree heights was measured.For that reason, height histograms for live and dead tree were created using the acquired LiDAR.In Figures 5a and 6a, the parameter height is calculated after the Digital Terrain Model (DTM) is subtracted but it includes an offset of around 20 m that corresponds to noise recorded from the LiDAR system after the pulse reached the ground.To further clarify this, the data are interpreted in DASOS [45] that defines the boundaries of the voxelised space according to the positions of the first and last waveform samples of each pulse.Full-waveform LiDAR systems usually record a predefined number of waveform samples.If a pulse reaches the ground before the predefined space for values is filled, the systems keeps digitising noise that correspond to locations below the ground.This recorded noise lies inside the boundaries of the voxelised space.When the DTM is subtracted, the boundaries still contain that noise, which appears as a constant offset from the ground on the height metrics of DASOS.Figures 5a and 6a show that a high percentage of locations indicating trees within the field data are classified as ground.To overcome this, in Figures 5a and 6a the histograms of heights are created using a 3 × 3 pixels window (corresponds to 2.4 m × 2.4 m) at each tree location of the field data.The highest value within the 3 × 3 window is taken.It is worth mentioning that as stated above a eucalypt tree reaches up to 40 m and this is also shown in the histograms considering that anything below 20 m is considered ground or below the ground due to the noise recorded.Those histogram depict the tree height differences within the dataset and guided the definition of the sizes of the 3D-windows used to extract information and train/test the classifier.

Methodology
The proposed methodology is a pipeline of processing steps that leads towards the detection of the dead standing eucalypt trees.Table 1 gives an overview of proposed pipeline and highlights the differences between the new multi-scale 3D-windows approach and the single size 3D-windows approach [53].The following sub-sections explain each step of the processing pipeline with more details.For an intensive explanation of the single 3D-windows approach and the overlapping steps (e.g., weighted k-Nearest Neighbour (k-NN) algorithm, seed growth algorithm), please refer to this article [53].

Voxelisation of Waveform Samples
In this paper, the waveform samples are voxelised using a 3D regular grid of resolution 0.8 m.Please note that the resolution of the probabilistic field created by k−NN (Section 2.2.4) was reduced to 1 m due to rounding errors occurs in the C++ implementation.The voxelisation process can be summarised as follow: 1.
The space is divided into a 3D regular grid; each cube is named voxel 2.
Noise reduction by applying low level filtering-in other words ignore waveform sample with low intensity profile since they most probably contain noise 3.
Each waveform sample -whose position is calculated as explained in Section 2.1.2-isassociated with the voxel that it lies inside 4.
Each voxel takes as value the average intensity of the waveforms that are associated with it 5.
The result is a discrete density volune (can be interpreted as a 3D grayscale image) with accumulated intensities of the waveforms More information about voxelisation, including visual explanation, were published by Miltiadou et al. [45].It further worth highlighting that while calculating the position of the waveforms, the digital terrain model is taken into consideration and the height of the terrain is subtracted before each waveform sample is associated with a voxel.

Extraction of Features from 3D-Windows
Composite information from 3D-windows (e.g., distribution of non-empty voxels) was extracted for characterising dead trees (positive samples), live trees (negative samples) and testing data using DASOS, an open source software whose main developer is the corresponding author of this paper [45].For cross validating the results and consequently avoiding overfitting and underfitting of the classifier, the 50 field plots (including field plots used twice due to overlapping flightlines) were shuffled and randomly divided into four validation datasets (Val1, Val2, Val3 and Val4).Each validation datasets contained 25 field plots for training the classifier and 25 field plots for testing the results of the proposed methodology.For the multi-scale 3D-windows approach, for each validation dataset we created three sets of positive, negative and testing samples; one for each window-size (Figure 7).Equal number of positive and negative samples were used to train the classifier to reduce bias.Dead trees were significantly less in number than live trees.A sample of live trees were carefully selected such that the distribution of heights was the same between dead and live trees.The validation dataset Val1 was used for defining the appropriate thresholds and filters of the methodology.
In the multi-scale 3D-windows approach the height categories are defined by observing the histograms in Figures 5 and 6.The trees listed within the training datasets are divided into three size categories: (1) CHM < 28, (2) 28 <= CHM < 38 and (3) CHM >= 38 (including the 20 m offset for the noise).The shape of the windows is a cylinder defined by its height and radius.Height and radius are divisible by 0.8, which is the length of a voxel.Using several tests, as reference [53] and the fact that the diameter should contain an odd number of voxel lengths for DASOS to work [45], the sizes of the windows are defined as follow: (1) height: 4.8 m, radius: 4.0 m, (2) height: 9.6 m, radius: 7.2 m and (3) height: 14.4 m, radius: 10.4 m.Nevertheless, these sizes could be modified.

Random Forest
The open source software DASOS, developed mainly by the first author of the paper, exported 28 features (e.g., distribution of non-empty voxels) for each training/testing sample.The k-NN algorithm is ideal for detecting objects with irregular shapes, since it does not fit them into a single statistical model.Nevertheless, it is sensitive to noise and, consequently, the random forest was used to reduce the dimensionality of the classifier from 28 to 10. Random forest generated multiple decision trees to identify the most important features [57], which were used in the k-NN algorithm.The ten most important features identified were: "Height_Std", "Top_Patch_Len_Std", "Dis_Std", "Per_Int_Above_Iso", "Top_Patch_Len_Mean", "Top_Patch_Len_Median", "Dis_Mean", "Dis_Median", "Sum_Int_Diff_Z", "Sum_Int_Diff_Y".Explanations about the extracted features are given at [45].

Weighted k-Nearest Neighbours Algorithm
The probability of each column-of the voxelised space-to belong to the dead tree population was calculated using the k-NN algorithm.Information about the area around a column of interest was extracted by placing a 3D-window with the column of interest to be in the middle; starting from the top of the voxelised space the window was moving downwards until the first non-empty voxel of the middle column was reached and the features are extracted from the area that lies inside the 3D window (Figure 8).For each column within the voxelised space, k-NN calculated the probability of this column to contain a dead tree and, therefore, reduced dimensionality from 3D to 2D; each pixel in the output 2D image is the possibility of a corresponding column to belong to the dead tree population.This 2D image is called probabilistic field [53].A feature vector is a vector containing many features extracted from a single 3D-window.Positive feature vectors were derived by extracting information from 3D-windows placed around dead trees and negative training samples were feature vectors derive from 3D-windows placed around live trees.It is worth highlighting that the same number of positive and negative training samples were extracted and used to train the classifier to reduce bias.Otherwise, a class confidence weighted k-NN algorithm could have been applied [58] to balance the high difference between the number of dead and live trees within the field data.For the single size 3D-windows approach one probabilistic field was created for each field plot, while for the multi-scale 3D-windows approach three probabilistic fields were created, one for each window size, and they were merged afterwards during the fusion step Section 2.2.5.Starting from the top of the voxelised space the 3D-window is moved downwards until the first non-empty column of the middle column is reached.Information are extracted from the area that lies inside the 3D-window (middle column and its surroundings).

Fusion of Probabilistic Field
The output of the k-NN algorithm was a probabilistic field, where each pixel contained the probability of a dead tree to be located at the geolocation of the pixel.In the single size 3D-windows approach only one probabilistic field was created, while in the multi-scale 3D-windows approach three probabilistic fields were created; one for each window size.The fusion step only applied for the multi-scale 3D-windows approach and it fused the three probabilistic fields into one.A script, developed for this study, took as input the three probabilistic fields and the CHM (with the 20 m offset).As shown in Figure 9, the CHM of each pixel defines from which of the three k-NN outputs the corresponding pixel value is assigned to the fused image.For example, if the CHM of pixel (x, y) is less than 28 m then the value of the fused image at pixel (x, y) will be equal to the value of the pixel (x, y) of the probabilistic field created using the small size window.

Location detection
The fused image for the multi-scale approach or the results of the k-NN for the single size 3D-windows approach was a probabilistic field showing the probability of each pixel to be dead.For noise reduction, median and averaging filtering was applied.Afterwards, pixels classified as ground (CHM < 20 m including noise recorded after the laser pulse hit the ground) and pixels with probability of been dead less than 60% were thresholded.The result was a set of segments, which were also filtered with median filtering.Examples of these steps are shown in Figure 10.
Finally, a seed growth algorithm was applied for identifying unique segments and assigning predicted locations of dead trees.The predicted locations were equal to the average pixel geo-location of each segment.The training data were divided into three categories according to their height (small: CHM < 28 m, medium: 28 m< CHM < 38 m and large CHM > 30 m).The k-NN algorithm is applied three times using a different 3D-window size each time.The maps in the middle show the three probabilistic fields created by k-NN: each pixel contains the probability to belong to the dead tree population.The fusion is done pixelwise.For instance, if the canopy height at pixel (x, y) belongs to the small population, then the pixel (x, y) of the final image takes the value of the probabilistic field created using the small size of 3D-windows and the trees with CHM < 28 m for training the classifier.

Cross Validation, Precision and Recall
A cross-validation approach was used to assess the performance of the proposed multi-scale 3D-windows methodology for tackling height variations and the results were compared with the single size 3D-windows approach.There were 50 field plots (unique associations of field data and flightlines) available for evaluation.As explained in Section 2.2.2, the field data were randomly divided four times creating this way four validation datasets (i.e., Val1, Val2, Val3 and Val4).Each validation dataset contained 25 field plots used for extracting positive and negative training samples.The other 25 field plots were used for testing the performance of the proposed methodology.
For correctly evaluating the accuracy of the results, the precision (i.e., how many of the locations that the classifier labelled as dead trees were correct) and recall (i.e., percentage of dead trees detected) were calculated.At first the following are counted:

•
True Positives (TP): how many dead trees the system detected correctly • False Positive (FP): how many locations the system wrongly labelled as dead trees • False Negatives (FN): how many dead trees were not detected Consequently • TP + FP: is the number of locations that the classifier predicted that a dead tree exist • TP + FN: is the number of dead trees that existed within the field data Precision is equal to (TP/(TP + FP)), while recall is equal to (TP/(TP + FN)).Precision and recall are negatively correlated; when thresholds are modified one is improved and the other one is worsened.Therefore, in order to prove that an algorithm works better than another one, it is important to calculate both and show that both of them have been improved.

Results
The single size 3D-windows and the multi-scale 3D-windows approaches were compared by calculating the precision and recall (Section 2.2.7) of each validation dataset.For evaluating the results of the prediction, a distance parameter was used and a predicted tree was considered as TP if an actual dead tree existed within that distance.The distance parameter was considered essential for the following reasons (1) there was unknown error within the field data and (2) a tree lay in an area, while the field data contained only locations of trees and the system predicted locations as well.Therefore, a predicted location may lie within an area that a dead tree covers, but not be equal to the location of a tree listed as dead in the field data.The distance parameter is, therefore, important.Considering the sizes of the windows and the sizes of the trees, a distance above 7 m is too far from the target.Additionally, even though 0.5 m was considered in the evaluation, it is a very short distance considering that the size of the pixels in the k-NN results were 1 m.This study investigates how the precision and recall are modified when the distance parameter is modified.Alternatively, a constant distance could have been chosen and this will have helped into investigating how precision and recall is modified by changing the filters applied.Tables 2 and 3 and Figure 11 contain the results of the single size 3D-windows approach.Tables 4 and 5 and Figure 12 show the results of the multi-scale 3D-windows approach.Figure 13 shows the average recall and precision of each implemented approach.On average, the precision and recall of the new methodology was improved by 2.03% and 27.6%, respectively.

Discussion
Dead wood is extremely important for managing biodiversity, since while decaying it provides resources for numerous organisms [1].Fungi plays a substantial role in the formation of hollows and wood decaying, which further supports biodiversity and, consequently, a resilient ecosystem [59].In Australia, tree hollows are formed by fungi and provide shelters to native arboreal mammals and birds [2,11].They are, therefore, important for managing biodiversity.For detecting snags from infrared imagery, Polewski et al. [60] used a two stages detections approach.At first Gaussian analysis was used to estimate locations of dead trees and they used prior knowledge about shapes and density of local areas manually labelled as snags.Similarly, Miltiadou et al. [53] used prior knowledge by extracting features from 3D-windows around dead and live trees, labelled in field data, and used these features to perform the detection.
The study presented in this article aims to increase resilience to tree heights and sizes, while working with native forests.It introduces the usage of multi-scale 3D-windows for detecting snags.For assessing the performance of the new proposed algorithm, the multi-scale 3D-windows approach has been compared with a single size 3D-windows approach and it was proven that it has improved detection of dead standing eucalypt trees without tree delineation.It was shown before that the single size 3D-windows approach performs better than a random distribution of predicted dead tree locations with equal average density per plot [53].Additionally, since the multi-scale 3D-windows approach improves both recall and precision once compared with the single size 3D-windows approach (Figure 13), the improvement in detecting dead standing eucalypt trees is doubtless (precision and recall are explained in Section 2.2.7).
Overall, the proposed methodology confers better results than what has been done before without tree delineation, but needs further improvements.It is noticeable that the recall of Val4 is lower when the multi-scale 3D-windows approach is applied, but its precision increases by 25.36%.This indicates potential over-prediction of dead trees; if too many dead trees are predicted, then precision is high since the algorithm manages to predict a high percentage of the actual dead trees but recall is low because the probability of predicted dead trees to be actual dead trees is low.This issue may be solved by adjusting the filters and thresholds.
Another limitation is the unknown noise within the field data and the limited number of trees measures in the field.The locations of the trees seem to have no consistency once plotted on CHM (Figure 2).This may occur because some of the dead trees are really small or sometimes only the trunk has remained, which is not acquirable by the LiDAR system.These data only produce noise and reduce the accuracy (precision and recall) of the classifier.This is a major limitation, since the classifier is trained using noisy field data, while at the same time there were occasions were less than 15 dead trees were used in the k-NN for a specific window size.Increased training data and manual inspection of them should improve prediction.
Despite the aforementioned limitations, it is demonstrated that on average the proposed multi-scale 3D-windows approach improves the prediction.This opens up possibilities of further research in native forests that contain complex shapes.The effectiveness of the proposed methodology could be tested in various applications, like biomass estimation [61], leaf area index [62] and estimating total stem volume [63].The sizes of the windows will have to be adjusted to the new site but this can be done by observing the distribution of tree heights using histogram analysis.
In this study, the k-NN algorithm was used because it does not fit the data into a single statistical model and it is, therefore, expected to give good results in detecting objects with variant shapes like snags.Nevertheless, Wu and Zhang [64] has recently showed that Support Vector Machines (SVM) work better on classifying tree species in relation to k-NN and this was occur due to the noise.This is reasonable considering that k-NN is very sensitive to noise.In the study presented here, the selected window sizes are smaller than the average height/size of each tree size category to reduce noise.Additionally, the random forest is first used to select the most important parameters for distinguishing dead from live trees.Only the most important parameters are used in the classifier for reducing dimensionality and, consequently, noise.Even though the selection of this approach is done based on knowledge about the acquired data, in future work it worth checking the performance of further state-of-art machine learning approaches.For example, SVM supports high dimensionality while classifying data and Zhao et al. [65] showed that SVM performs better than the maximum likelihood classifier and linear regression models for estimating various forest parameters.Furthermore, with the increased computational capabilities of computers neural networks are used in recent literature for advancing forest inventories [66,67].Working with small classification datasets there is a risks that neural networks may conclude that the relevant features within the training dataset are noise.Nevertheless, this possibility has been decreased with the big data era and increased training samples.
It further worth highlighting the important structural parameters identified for distinguishing dead trees from alive, which are (1) Height_Std: the standard deviation of the heights inside the 3D-window, (2) Top_Patch_Len_Std the standard deviation of the length of all the top patches (a top patch is defined as the number of adjacent non-empty voxels starting counting from the top of the column of interest) and (3) the Standard Deviation of the distances between the central voxel and every non-empty voxel lying inside the 3D-window.These parameters are reasonable considering that dead trees have less leaves in comparison to live trees and are, therefore, more likely to have bigger height differences within the 3D-windows.These parameters can be adopted in other sites and applications; for example for the detection of infected trees and for distinguishing deciduous from evergreen trees during the leaf-off season due to their structural similarities.

Conclusions
This study introduces the usage of multi-scale 3D-windows for tackling height variations while detecting dead standing eucalypt trees, which are important for managing biodiversity in native Australian forests.In the proposed methodology, the LiDAR waveform samples are voxelised by inserting them into a regular 3D-grid.Afterwards, by observing the height histogram of eucalypt trees, the trees are divided into three height categories.For each category, an appropriate size of 3D-window is used for extracting features characterising dead (positive samples), live (negative samples) trees and columns of the voxelised space (testing samples).This information is used for training and testing a multi-scale 3D-windows classifier.The results of the proposed methodology were cross-validated and compared with the results of a single size 3D-windows approach.The precision and recall of the new methodology was improved by 2.03% and 27.6%, respectively.As explained in Section 2.2.7, this proves that the prediction is improved with the multi-scale approach, since both precision and recall are increased.Consequently, the new methodology can be used for tackling height variations while detecting dead standing Eucalypt trees in native Australian forests.Furthermore, this study opens ups possibilities of new research directions and applications related to the proposed methodology for deriving forest related parameters (e.g., crown detection, biomass estimation and leaf area index) from native forests with variant tree heights and sizes.

Figure 1 .
Figure 1.The study area used is depicted by red.The entire study area was scanned.The field plots used for evaluation were randomly scattered within that area and their locations are depicted with green circles.

Figure 2 .
Figure 2. Example of a field plot for understanding the complexity of detecting dead eucalypt trees.The image on the left depicts the locations of both live (red stars) and dead (yellow stars) trees on Canopy Height Model (CHM).The image on the right shows the airborne laser scanning data, created by processing waveform data, and the estimated locations of the dead trees are shown in circles.

Figure 3 .
Figure 3. Example of acquired LiDAR data (discrete returns extracted from waveform data) depicting how the system failed to acquire information about tree trunks [53].

Figure 4 .
Figure 4. Illustration of how full-waveform LiDAR systems work.They emit a laser pulse and digitise the entire backscattered signal.Digitisation is done at equal time space intervals and the average distance between two coherent waveform samples is constant.
(a) Live trees (b) Dead trees Figure 5. Histograms of elevation at tree locations listed within the field data using a single column to measure the elevation.Histogram (a) includes the heights of the Live Trees and histogram (b) the heights of the dead trees.
(a) Live trees (b) Dead trees Figure 6.Histogram of the highest elevation within 3 × 3 windows of columns with centre the location of trees listed within the field data.Histogram (a) includes the heights of the Live Trees and histogram (b) the heights of the dead trees.

Figure 7 .
Figure 7. Feature vectors extracted for training and validating the multi-scale 3D-windows methodology.

Figure 8 .
Figure 8. Features from local areas are extracted using 3D cylindrical windows.For each column of interest, a window is placed such that the column of interest lies in the middle of the 3D-window.Starting from the top of the voxelised space the 3D-window is moved downwards until the first non-empty column of the middle column is reached.Information are extracted from the area that lies inside the 3D-window (middle column and its surroundings).

Figure 9 .
Figure 9. Pixelwise fusion of the three probabilistic fields into one using the CHM as a reference.The training data were divided into three categories according to their height (small: CHM < 28 m, medium: 28 m< CHM < 38 m and large CHM > 30 m).The k-NN algorithm is applied three times using a different 3D-window size each time.The maps in the middle show the three probabilistic fields created by k-NN: each pixel contains the probability to belong to the dead tree population.The fusion is done pixelwise.For instance, if the canopy height at pixel (x, y) belongs to the small population, then the pixel (x, y) of the final image takes the value of the probabilistic field created using the small size of 3D-windows and the trees with CHM < 28 m for training the classifier.

Figure 10 .
Figure 10.Processing pipeline: (a) contains the results of the k-NN algorithm or the fused image, (b) the result of applying Median filtering on the previous image, (c) result of applying an averaging filter on the previous image, (d) result once the ground is thresholded, (e) results once the k-NN threshold is applied (f) Median filtering.

Figure 11 .
Results of the algorithm using a single size window: (a) Recall, (b) Precision.

Figure 12 .
Results of the algorithm using three windows with different sizes: (a) Recall, (b) Precision.

Figure 13 .
Average results of the two approaches: (a) Recall, (b) Precision.

Table 1 .
List of processing step with associated sections.

Table 2 .
Recall while using a single size of 3D-window.

Table 3 .
Precision while using a single size of 3D-window.

Table 4 .
Recall while using multi-scale 3D windows.

Table 5 .
Precision while using multi-scale 3D windows.