A graph cut approach to 3D tree delineation, using integrated airborne LiDAR and hyperspectral imagery

2.8 MB·

by  Juheon Lee

A graph cut approach to 3D tree delineation, using integrated airborne LiDAR
and hyperspectral imagery
JuheonLeea,b,DavidCoomesa∗,Carola-BibianeScho¨nliebbXiaohaoCaia,b,JanLellmannb,MicheleDalpontea,c,
YadvinderMalhid,NathalieButtd,e,MikeMorecroftf,
aForestEcologyandConservationGroup,DepartmentofPlantSciences,
UniversityofCambridge,CB23EA,UK
bImageAnalysisGroup,DepartmentofAppliedMathematicsandTheoreticalPhysics(DAMTP),
UniversityofCambridge,CB30WA,UK
cDepartmentofSustainableAgro-ecosystemsandBioresources,ResearchandInnovationCentre,FondazioneE.Mach,ViaE.Mach1,38010San
Micheleall’Adige(TN),Italy
dEnvironmentalChangeInstitute,SchoolofGeographyandtheEnvironment,
UniversityofOxford,OX13QY,UK
7 eCentreforBiodiversityandConservationScience,TheUniversityofQueensland,
1 StLucia,4072,Qld,Australia
0 fNaturalEngland,CromwellHouse,15AndoverRoad,Winchester,SO237BT,UK
2

n
a
J

4
2 Abstract

]
V Recognisingindividualtreeswithinremotelysensedimageryhasimportantapplicationsinforestecologyandman-
C
agement. Several algorithms for tree delineation have been suggested, mostly based on locating local maxima or
.
s
c invertedbasinsinrastercanopyheightmodels(CHMs)derivedfromLightDetectionAndRanging(LiDAR)dataor
[
photographs. However, these algorithms often lead to inaccurate estimates of forest stand characteristics due to the

1
v limited information content of raster CHMs. Here we develop a 3D tree delineation method which uses graph cut
5
1 todelineatetreesfromthefull3DLiDARpointcloud, andalsomakesuseofanyopticalimageryavailable(hyper-
7
6 spectralimageryinourcase). First,conventionalmethodsareusedtolocatelocalmaximaintheCHMandgenerate
0
. aninitialmapoftrees. Second,agraphisbuiltfromtheLiDARpointcloud,fusedwiththehyperspectraldata. For
1
0
computational efficiency, the feature space of hyperspectral imagery is reduced using robust PCA. Third, a multi-
7
1
class normalised cut is applied to the graph, using the initial map of trees to constrain the number of clusters and
:
v
i theirlocations. Finally,recursivenormalisedcutisusedtosubdivide,ifnecessary,eachoftheclustersidentifiedby
X
r the initial analysis. We call this approach Multiclass Cut followed by Recursive Cut (MCRC). The effectiveness of
a
MCRCwastestedusingthreedatasets: i)NewFor, whichincludesseveralsitesintheAlpsandwasestablishedfor
comparingsegmentationalgorithms, ii)aconiferousforestintheItalianAlps, andiii)adeciduouswoodlandinthe
UK.TheperformanceofMCRCwasusuallysuperiortothatofotherdelineationmethods,andwasfurtherimproved
byincludinghigh-resolutionopticalimagery. SinceMCRCdelineatestheentireLiDARpointcloudin3D,itallows
individualcrowncharacteristicstobemeasured. MCRCiscomputationallydemandingand,likecurrentCHM-based
approaches,isunabletodetectunderstorytrees. Nevertheless,bymakingfulluseofthedataavailable,graphcuthas
1
thepotentialtoconsiderablyimprovetheaccuracyoftreedelineation.
Keywords: Treesegmentation,remotesensing,LiDAR,hyperspectralimage,opticalimagery,normalisedCut.
1. Introduction
There is much interest in using remote sensing to map individual tree crowns (ITCs) (Heinzel and Koch, 2012;
Dalponte et al., 2014) and measure various attributes of the identified trees. Until recently, applications of remote
sensingdatatovegetationmonitoringhavefocussedmostlyonproducingrasterised2Dmaps,witheachpixelsum-
marizing information from the many individual plants within them (Clark et al., 2005; Asner and Martin, 2008).
For example, NASA’s Landsat 8 satellite produces imagery at 30m spatial resolution, which is too coarse to detect
individualtrees,buthasbeenusedtomapforesttypes,deforestation,blahblahamongothers[[Hansen]]andpixel-
basedapproachestoanalysingairborneremotesensingarealreadyintegratedintoanationalforestinventoryprogram
(e.g. Finland (Tomppo, 1993)). However, tree-centric approaches have the potential to advance forest research, by
keeping track of individual responses to pest and pathogen outbreaks, selective logging, fire, invasive species and
climate change (Asner et al., 2008b; Andersen et al., 2014; van Ewijk et al., 2014). Earth observation technology
is producing information at increasingly high spatial resolution, making ITC approaches an attractive alternative to
pixel-basedmethods. Inparticular,airborneLiDARproducesa3Dpointcloudindicatingwherelaserpulsesemitted
from the transceiver have reflected off leaves, branches and the forest ground, making it possible to map individual
treesovertensofthousandsofhectares(Chenetal.,2006;HeinzelandKoch,2012;Colganetal.,2012). Inaddition,
airbornehyperspectralimagerycanbeusedtoestimatethephysicalandchemicalpropertiesofcanopies, andwhen
usedalongsideLiDARcanmapthesepropertiesatanindividualtreelevel(Asneretal.,2007). This”spectranomic”
approachhasbeenusedtomapinvasivetreespeciesinHawaiianrainforests(Asneretal.,2008a,b)andtoquantify
thespatialvariationinbiochemicaldiversityintropicalregions(AsnerandMartin,2008,2011).
Tree-centric approaches are not widely adopted by the forest remote sensing community as yet, in part because
of problems with accuracy. Classical delineation approaches work with rasterized canopy height models (CHMs)
derived from the LiDAR point cloud. These methods include watershed algorithms (Chen et al., 2006; Koch et al.,
2006; Kwak et al., 2007; Yu et al., 2011), variable window filtering (Hyyppa et al., 2001; Solberg et al., 2006),
multi-scaleedgesegmentation(Brandtbergetal.,2003), andattentivevisionmethods(Palenichkaetal.,2013). All
these approaches share several problems: (a) the smoothing process determines the number of trees detected in the
CHM:toostrongasmoothingfactorleadstounder-segmentation, whiletooweakafactorgeneratesmanytree-like
∗DavidCoomes,[email protected]
PreprintsubmittedtoElsevier January25,2017
artifacts(Maltamoetal.,2014);(b)sub-canopytreesareimpossibletodetectastheyallrelysolelyoncanopysurface
geometry;and(c)theinterpolationandsmoothingprocessesinvolvedingeneratingCHMsresultinunderestimation
of tree heights, meaning that additional post-processing is needed to rectify the results (Solberg et al., 2006; Koch
etal.,2006). Inordertoaddresstheseproblems,moreadvancedmethodsthatexploittheentireLiDARpointcloud
have been developed. These methods include k-mean clustering (Morsdorf et al., 2004), normalised cut (NC) (Shi
andMalik,2000;VonLuxburg,2007;Reitbergeretal.,2009;Yaoetal.,2012),adaptiveclustering(Leeetal.,2010),
support vector machine (SVM) (Secord and Zakhor, 2007; Zhao et al., 2011), and exploiting the spacing between
top of trees (Li et al., 2012). Most of these approaches were developed for managed coniferous forests, which are
relatively straightforward to delineation because conical crowns have well defined peaks and forest size structure is
simple. Benchmark datasets available to compare approaches also focus on coniferous forests (NEWFOR, 2012).
Much less work has been done in tropical or temperate broadleaf forests, where intermingled dome-shaped crowns
makedelineationmorechallenging(Reitbergeretal.,2009;Yaoetal.,2012;Lietal.,2012;NEWFOR,2012;Heinzel
andKoch,2012;Immitzeretal.,2012;Colganetal.,2012). TheapproachofDuncansonetal.(2014)holdpromisein
thisregard;theyfirstdelineatedtreesintheuppercanopyusingawatershedapproach,thensearchedfor”troughs”in
theverticalstructureofthelocal3Dpointcloudthatallowedthemtostripawaythetallertreesandusethewatershed
algorithm a second time, to delineate subcanopy trees (Duncanson et al., 2014). In principle, the fusion of high
resolutionopticalimagerywithLiDARdatashouldleadtoimprovementsinITCdelineation(Chenetal.,2006;Koch
et al., 2006; Kwak et al., 2007; Hyyppa et al., 2001; Solberg et al., 2006; Brandtberg et al., 2003; Dalponte et al.,
2011;Yuetal.,2011)byhelpingtodistinguishneighbouringtreesthroughdifferencesintheirradiometricproperties
(Kaasalainenetal.,2009;Korpelaetal.,2010b,a). Aerialphotographsandmulti/hyperspectralimagerycouldallbe
usedforthispurpose,aslongastheirspatialresolutionishighenough(i.e. thepixelsizeissmallerthantheminimum
crownsizethatweneedtodetect)(Koch,2010;Sua´rezetal.,2005;Holmgrenetal.,2008;Breidenbachetal.,2010;
Colgan et al., 2012; Heinzel and Koch, 2012; Dinuls et al., 2012; Jakubowski et al., 2013). However, multi-sensor
approaches are only possible if the different data are accurately co-aligned, thus image registration must be applied
prior to their fusion (see (Dawn et al., 2010; Le Moigne et al., 2011)). A second issue is that extracting feature
informationdirectlyfromhighdimensionaldata-suchasthehyperspectraldatasets-oftenleadstoinaccurateresults
(Dalponte et al., 2008). Therefore, dimensional reduction is required before applying any delineation algorithm,
usingfeatureextractiontechniquessuchasprincipalcomponentanalysis(PCA)(Cande`setal.,2011),orbyselecting
influentialfeaturesfromtheoriginalbands[[check]](Dalponteetal.,2008,2011).
ThisstudyseekstoovercomesomeoftheissuesassociatedwithITCdelineation,bydevelopinganewapproach
graphcutapproach,basedonNormalisedCut(ShiandMalik,2000;VonLuxburg,2007;Reitbergeretal.,2009;Yao
3
etal.,2012),whichcombinesopticalandLiDARdata. NormalisedCutisawell-establishedapproachforgrouping
pointsand/orpixelsintodisjointclusters. Itstartswithamatrixofsimilaritymeasuresbetweenallpossiblepairsof
pointsand/orpixels,andusestheeigenvectorsofthatmatrixtodistinguishgroups(ShiandMalik,2000). Inthecase
of LiDAR data, the similarity matrix is derived from the physical distance between points (nodes) in 3D space. In
thecaseofhyperspectraldata,thematrixisderivedfromtheirradiometricsimilarityandphysicaldistancesbetween
pixels(nodes)in2Dspace. NCseekstopartitionthegraphintoclusterswithhighsimilaritiesbetweenthenodesof
thesameclustersandalowsimilaritiesbetweennodesfromdifferentclusters. TheadvantageoftheNormalisedCut
approachisthatgraphweightscanbedefinedusingopticalimageryalongsideLiDAR,thusprovidingaframework
forfusingdifferenttypesofremotesensingdatasets.
Ourmainobjectiveistodescribe,andevaluate,agraphcutapproachthatcanbeusedtodelineateITCsdirectly
froma3DLiDARpointcloud, usingsupplementaryinformationforopticalsensors. Todothiswe(a)describethe
data processing pipeline, including an efficient way of fusing LiDAR data and optical imagery and various graph
cut methods ; (b) examine the capability of MCRC to detect understory trees and correctly segment canopy trees
by working with forest plot data from coniferous and broadleaf woodlands. The paper is organized as follows: in
Section2,thegeneralmathematicalprinciplesofthenormalisedcutapproachareoutlined. Theapplicationofthese
principlestotreedelineationinwoodlandsisintroducedinSection3.Thetestdatasetsusedtoexemplifyourapproach
aredescribedinSection4. TheperformanceofourapproachisevaluatedinSection5. Section6discussesMCRCin
relationtootherapproaches,andgivesrecommendationsforfuturework.
2. Thegeneralprinciplesofthenormalisedcutapproach
This section provides a formal outline of the normalised cut approach (Shi and Malik, 2000; Reitberger et al.,
2009). A graph G is a pair of sets, G = (V,(cid:15)), where V is the set of N vertices and (cid:15) is the set of edges. Each
edgew ∈ (cid:15) correspondstoanon-negativesimilarityweightbetweentwoverticesi, j ∈ V. Theobjectiveofbinary
ij
graphcutistopartitionthegraphintotwodisjointsetsAand Bbycuttingedgesthatconnectthetwosets,suchthat
A∪B=V andA∩B=∅. WedefinethecutasthesumovertheweightsofalledgesthatconnectAandB,thatis
(cid:88)
cut(A,B)= w (1)
i,j
i∈A,j∈B
4
(cid:80)
Wedefineassoc(A,V)tobethetotalweightofconnectionsfromnodesinAtoallnodesinthegraph(i.e., w ).
i∈A,j∈V ij
ThenormalisedcutmethodfindssetsAandBbyminimisingthefollowingenergyterm:
cut(A,B) cut(A,B)
Ncut(A,B)= + . (2)
assoc(A,V) assoc(B,V)
Inordertosolve(2),tofindsolutionx∈RN ,wereformulateitas:
minxTD−12(D−W)D−12x
x∈RN (3)
s.t. xTD121=0, xTx=1,
where D ∈ RN×N is a diagonal matrix with diagonal entries d = (cid:80)N w , W ∈ RN×N is a symmetric matrix with
i j=1 ij
entitieswij,and1∈RN isanall-onesvector. Solutionsarefoundbycalculatingtheeigenvectorsofmatrix D−21(D−
W)D−21 . The smallest eigenvector is D211 , which is a trivial solution, and is ignored. It is the second smallest
eigenvectorthatistakenasthesolution. V isthensplitintotwosetsbythresholdingthesolution x,forexampleby
takingthemeanvalueofx.
Theapproachcanbeextendedtosearchformultipleclasses,byapplyingthebinarygraphcutrecursively(i.e.sets
AandBarefurthersubdividedintofoursets,andeachofthemmaybesubdivided,andfurthersubdivided)untilthe
processisterminatedbyastoppingrule.ThedecisionastowhetherornottomakeasplitdependonwhethertheNcut
energy value exceeds some predetermined threshold (Shi and Malik, 2000). By this recursive application of graph
cut,theindividualtreesintheforestcanbedelineated. However,therecursiveschemeiscomputationallyinefficient
becauseitneedstosolveequation(3)repeatedlyuntilitreachesthispredefinedthreshold. SinceLiDARdatacontains
millionsofpointsperhectare,recursivegraphcutrequireshugecomputationalpowertoworkondatasetslargerthan
a few square metres. A further issue with the recursive approach is that equation (3) uses only the second smallest
eigenvector (Shi and Malik, 2000; Von Luxburg, 2007), discarding information from subsequent eigenvectors that
could help refine the partitioning. Finally, it is difficult to incorporate priors (i.e. initial guesses of the location of
clusters) using this approach, which turns out to be important when delineating trees (see later). For these reasons,
there are advantages to using a multiclass normalised cut approach, that searches for a predetermined number of
classes, instead of using recursive binary cut (Von Luxburg, 2007). The multiclass problem can be understood as
follows: let the solution matrix X = (x ,··· ,x ) ∈ RN×C whereC be the number of clusters. Then, the multiclass
1 C
5
problemcanbeexpressedinasimilarwaytoproblem(3):
min tr(cid:0)XTD−12(D−W)D−21X(cid:1)
X∈RN×C (4)
s.t. xTi D121=0, xTi xi =1, i=1,...,C,
whereV issplitintoC setsbyeitherk-meansorspectralrotation. Thisapproachiscomputationallyefficient, since
thenumberofclustersisfixedatC,equation(4)needstobesolvedonlyonce. However,multiclassnormalisedcut
hastwoproblemswhenappliedtoforests.Thefirstproblemisthatthenumberoftreeshastobesetinadvance,which
somewhatdefeatsthepurposeoftreedelineation! Thisissueisresolvedbytakingeachoftheclustersidentifiedby
MCandapplyingabinarynormalisedcuttothemrecursivelytoidentifyfurthertreeswithinthecluster. Thesecond
issue is that preliminary trials showed that MC performs poorly unless the algorithm is given some clues as to the
whereabouts of trees. To resolve this problems, we first estimate the locations of tree tops from the local maxima
of the CHM and use these locations as priors, providing method (4) with an estimate of the number of clusters and
theirpositions. Constrainednormalisedcuthasbeenproposedby(Huetal.,2013)buthasneverbeenusedforITCs
delineation. This scheme regards a prior as an additional constraint to the solution of (4), minimising cut energy
but also satisfying the condition that the correlation between the solution and the prior is larger than or equal to a
predefinedvalue(κ).
Formally, let S = (s ,··· ,s ) ∈ RN×C be a matrix of priors. Then the MultiClass Normalised Cut with Priors
1 C
(MC)approachisgivenby
min tr(cid:0)XTD−12(D−W)D−12X(cid:1)
X∈RN×C (5)
s.t. xTi D121=0, xTi xi =1, xTi si ≥κ, i=1,...,C,
where κ is a correlation parameter. The solution of equation (5) gives C separate clusters of data. The correlation
term is a hard constraint, which must be satisfied. In other words, the solution must have C non-empty disjoint
clusters. Thismethodismuchfasterandmoreefficientthansolvingbinaryclusteringrecursivelybecausethenumber
ofclustersisfixedandequation(5)issolvedjustonce.
3. Methods
The data processing pipeline shown in Figure 1 has six steps: A. LiDAR data is separated into ground returns
from which a digital elevation model (DEM) is constructed, and object returns, from which a canopy height model
(CHM) is constructed; B. if optical imagery is available, a state-of-the-art feature reduction method - robust PCA
6
Figure1: TheworkflowusedtodelineateindividualtreecrownfromLiDARdata(solidline),andfromLiDARdatafusedwithopticalimagery
(dashedline).
(rPCA)(Cande`setal.,2011)-isusedtoreducethenumberofhyperspectralfeatureswithintheco-aligneddataset,
to speed up processing; C. if optical imagery is available, LiDAR and optical imagery are registered (precisely co-
aligned)usingtheNGF-Curvmethodthatwedevelopedpreviously(Leeetal.,2015);D.aconventionaldelineation
approach based on the CHM is used to identify likely locations of upper-canopy trees; E. These locations are as
priorsamulticlassnormalisedcut(MC);F.EachoftheclustersrecognisedbyMCaresubjectedtorecursivebinary
cutting. WecallthisMCRC(MultiClassNormalisedCutwithPriorsfollowedbyRecursiveNormalisedCut). Note
thatthismethoddelineatesITCsdirectlyfromthe3DLiDARpointcloud,soITCsarenotinfluencedbyinterpolation
orsmoothingerrorsprevalentinCHM-basedapproaches. InthefollowingsectionweexplaineachstepinFigure1in
moredetail.
A. Ground vs object filtering of point cloud We performed initial modelling of terrain and canopy heights
from the liDAR datasets using Tiffs 8.0: Toolbox for LidAR Data Filtering and Forest Studies, which employs a
computationally efficient, 25 grid-based morphological filtering method described by Chen et al. (2007). Outputs
includedfilteredgroundandobjectpoints,aswellasdigitalterrainmodels(DTM)andcanopyheightmodels(CHM).
B.FeatureextractionuserPCAHyperspectralimageryisinformationrich-oneofourdatasetshasinformation
collected from 361 contiguous wavebands. Using such a highly dimensional data in a graph cut is computationally
intensive, and making it practically difficult to exploit the information (Dalponte et al., 2008). To alleviate this
problem, the rPCA feature reduction technique is used in order to reduce the high dimensionality features space
7
toafewmeaningfulfeatures.ConventionalPCAissensitivetonoiseindata.Incontrast,rPCAisdesignedtorobustly
recoveralowrankmatrixLfromacorruptedmeasurementmatrix M,toleaveasparsematrixofoutliersS (Cande`s
etal.,2011). rPCAcanberepresentedasthefollowingminimizationproblem:
min{rank(L)+λ(cid:107)S(cid:107) } s.t. M = L+S,
0
L,S
where(cid:107)·(cid:107) isthel -normwhichimposesasparsitypropertyonS,rank(·)isthedimensionsofvectorspacesspanned
0 0
bythecolumnsorrowsofamatrix,andλisaregularisationparameter. Sincethisoptimisationproblemisintractable,
ingeneral,therankandthel -normareusuallyreplacedbythenuclearnorm(cid:107)·(cid:107) (sumofsingularvalues)andthe
0 ∗
l -norm(sumoftheabsolutevaluesofthewholeentries)respectively. Thisresultsinthefollowing:
1
min{(cid:107)L(cid:107) +λ(cid:107)S(cid:107) } s.t. M = L+S. (6)
∗ 1
L,S
This objective function is convex, so it can be solved by various convex optimisation algorithms. In this paper,
thealternatingdirectionmethodofmultiplierswasused(YuanandYang,2009). Weextractedthelowrankparts L
whichcorrespondstotheprincipalcomponentsinclassicPCA.Thefirstprincipalcomponentwasignoredbecauseit
containedilluminationinformationratherofusefulfeaturesofITCs(Tochonetal.,2015).Thesecondtofifthprincipal
componentswereextractedandassignedtocorrespondingLiDARpointsbyusinghorizontalgeospatialcoordinates.
IfthereismorethanoneLiDARpointinapixelofhyperspectralimagery,thenallpointsinthepixelwereassigned
thesamerPCAcoefficient.
C. Registration of remote sensing datasets LiDAR data and hyperspectral imagery are not usually precisely
co-alignedwhendeliveredbythedataprovider. Cameradirection,topographyandlensdistortionallaffectthequal-
ity of hyperspectral imagery, and LiDAR boresight is usually more accurate than that of the hyperspectral sensor,
so inaccuracies remain even after geometric correction. To co-align these data, registration of LiDAR and optical
imagerycanbeconductedusingNGF-Curvalgorithm,asproposedin(Leeetal.,2015). Thisnon-parametricregis-
trationmethodusesnormalisedgradientfieldsimilaritymeasureswithcurvatureregularisation. Comparedwiththe
traditionalparametricregistrationmethods(e.g.,(LeMoigneetal.,2011;Lietal.,2009)),theNGF-Curvmethodcan
handlenonlineardistortionandco-alignmulti-sensorimagerywithoutanygroundcontrolpoints. Thedetailsofthis
methodaredescribedin(Leeetal.,2015)andreferencestherein.
D.LocalmaximadetectionLocalmaximawithintheLiDARpointcloudprovidethepriorinformationontree
locationsinthispaper. ThoselocalmaximacanbeeasilyextractedfromtherasterizedCHM,usingamovingwindow
approach (Hyyppa et al., 2001) or a watershed approach (Chen et al., 2006). We used a marker-based watershed
8
approachfortreetopdetectionimplementedinTIFFS(Chenetal.,2006),comparingitefficacywithotherapproaches
usingtheNewForbenchmarkdataset(seebelow). AllLiDARpointswithin0.7mradiusofeachlocalmaximumwere
identifiedasbelongtothesamecluster. Weusedtheseclustersaspriors,thusenforcingthesolutionofequation(5).
Themarker-basedwatershedapproachisjustoneofthepossiblemethodstosetuppriors(Reitbergeretal.,2009)(see
Section5).
E.MultiClassNormalisedCutwithPriors(MC)TobuildthegraphforMC,weightsneedtobeassignedtothe
vertices,whicharegivenbytheLiDARpoints. WeusedanormalisedweightthatisafunctionofEuclideandistances
between vertices i and j in horizontal (x,y) and vertical (z) space, as well as the similarity of their hyperspectral
features(fts):
(cid:107)(xy)i−(x,y)j(cid:107) (cid:107)zi−zj(cid:107) (cid:107)ftsi−ftsj(cid:107)
wij = e σ2xy ×e σ2z ×e σ2fts , (7)
wherebandwidthparametersσ ,σ ,σ acttonormalisethefunction,andareparametersselectedbytheuser.
xy z fts
Forconstructingthegraph,weobservethatafullyconnectedgraphrequiresO(n2)memorycomplexity,whichisnot
practical. Instead,ad-neighbourhoodsamplingstrategyisadopted,whereweightsarecomputedonlywithinaradius
d of a vertex. In our examples, d ranged from 0.5m to 2m depending on the point density of LiDAR (lower radii
athigherdensitiestoreducethememorycosts). Equation(5)wassolvedwithad-neighbourhoodsimilaritymatrix
andpre-defined clusterstaken fromthe localmaxima. TheMC approachsegments the3D LiDARpoint cloudinto
thesamenumberoftreecrownsasidentifiedbytraditionalCHM-basedmethods,becausethisinformationisusedas
a”hard”prior. Italsosuffersfromthesameproblemsasclassicapproachesintermsoffailingtodetectunderstory
trees.
F. Recursive Normalised Cut (RC) The RC method (3) described in Section 2 is effective at ITC delineation,
including the detection of understory trees(Reitberger et al., 2009), but is computationally costly if applied to the
whole dataset. For this reason, we applied RC to each of the clusters obtained from an initial MC, to provide an
opportunityforcanopy”trees”tobefurthersubdividedandforsubcanopytreestobedetected.
4. DatasetDescriptionandDesignofExperiments
TheaccuracyoftheMCRCalgorithmwastestedon(a)asetforestplotslocatedintheAlpswhichformpartof
theNewForbenchmarkingproject,establishedspecificallyforthepurposeofcomparingITCalgorithms(NEWFOR,
2012),(b)aconiferousforestlocatednearTrentointheItalianAlps,and(c)alowlanddeciduousforestlocatednear
Oxford,UK.
(a) The NewFor LiDAR Single Tree Detection Benchmark Dataset consists of LiDAR and ground-truth in-
9
formationfrom14surveysitesintheAlps(10pilotareasin6countries)(Eysnetal.,2015). Amajoradvantageof
workingwiththeNewForbenchmarkdatasetisthatitprovidesanobjectivemeansofcomparingourapproachwith
others, and includes sophisticated validation software with which to evaluate algorithms by matching ITCs derived
fromLiDARwithknowntreelocationsinthefield. Thegroundtruthdatawereprovidedwithgeocoordinates, tree
height, DBH and canopy volume information. The errors of geocoordinates were less than 1m. The LiDAR point
densitywasmorethan10perm2 in12outof14studysite. TherangesoftheLiDARpointdensityintheNewFor
benchmarkdatasetwerefrom4to121perm2. AdisadvantageoftheNewFordatasetwithregardtoourproposedde-
lineationprocedureisthatitdoesnotcontainanyopticalimagery(i.e. weworkedwiththepipelineshownwithsolid
linesinFigure1). Notethatthesedatasetsareprimarilyconiferous,whicharerelativelystraightforwardtodelineate
becauseconifershavedistinctpeakstotheircrowns.
(b) The Italian Alps dataset was collected from a location near Trento. It consists of hyperspectral imagery,
LiDAR data and ground-based tree maps for 7 plots dominated by coniferous trees. Each plot is a circle of 15m
radius. Intheseplots,alltreeswithDBHabove1cmwereaccuratelygeoreferencedbydifferentialGPSandmanually
correctedwithlocalreferencetreesfromLiDARdata. Theestimatederrorofthegroundtruthoftreepositionswas
1m. The hyperspectral imagery were collected with an AISA Eagle sensor, covering 400–970nm with 61 spectral
bands,whiletheLiDARdatawereacquiredbyaRieglLMS-Q680isensoratanunusuallyhighpointdensity(≥87
pointsperm2). Hyperspectraldatawerecollectedon13th ofJune2013, whileLiDARdatawerecollectedbetween
7thand9thofSeptember2012.
(c) The English broadleaf woodland dataset was collected from Wytham Woods, Oxfordshire, England. It
containshyperspectralandLiDARdataover18hectaresoftemperatewoodlanddominatedbydeciduousangiosperm
species. The plot is fully mapped and all trees with DBH above 5cm are permanently tagged. As tree height was
measuredphysicallyonlyforsomeselectedsamplesweusedallometricequationstoestimatetreeheightforallthe
othertrees. Theestimatedpositioningerroroftheplotcornersisapproximately2m,whiletreepositionsarelocated
within about 5m. Hyperspectral imagery was collected in June, 2014 using an AISA Fenix sensor by the airborne
research and survey facility of the national environmental research council of UK (NERC-ARSF). It covers 400–
2500nm with 361 spectral bands. LiDAR data were collected by a Leica ALS-50 II scanner simultaneously with a
AISAFenixhyperspectralimagery. TheLiDARpointdensitywas6pointsperm2.
TheoptimalparametersfortheMCRCwerefoundbytrialanderror(Table1).
The validation of the ITC delineation was conducted using the tree matching software provided by the NewFor
project(Eysnetal.,2015;Kaartinenetal.,2012),whichcomparesrelativepositionsandheightsofsegmentedtrees
withthoserecordedinthegroundplots.Specifically,itmeasures2DEuclideandistanceandheightdifferencebetween
10