forked from edinburgh-university-OOSA/active_sensing
-
Notifications
You must be signed in to change notification settings - Fork 0
/
practical_week10.txt
135 lines (62 loc) · 6.78 KB
/
practical_week10.txt
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
Plain text instructions for week 10.
These can be copied and pasted in to the terminal to run the commands (word reformats the dashes to make them unrunable).
Install programs:
The command for the simple install is:
>> /geos/netdata/avtrain/data/3d/active_sensing/week10/active_sensing/installGedi.bash
>> source ~/.bashrc
To test the software has successfully compiled, type:
>> gediMetric -help
Viewing data
A script to plot the waveforms in the HDF5 file is available in:
/geos/netdata/avtrain/data/3d/active_sensing/week10/active_sensing/gediHandler.py
Do not worry about the contents of the script. It can be run from the command line without needing to modify the code (though you may modify if you want to).
Set up a shortcut to that directory:
>> bin="/geos/netdata/avtrain/data/3d/active_sensing/week10/active_sensing"
>> dataDir="/geos/netdata/avtrain/data/3d/active_sensing/week10/gedi_sonoma_area"
The options can be seen by typing:
>> python3 $bin/gediHandler.py --help
To use it to plot all waveforms in a given HDF5 file, type:
>> python3 $bin/gediHandler.py --input $dataDir/gediSonoma.power.1.day.200.0.h5 --outRoot $outRoot
Where “$filename” is the HDF5 file and “$outRoot” is a word to label the output files as (eg. “test1”). It will produce png files of the waveforms with the ground overlayed. Look at a few waveforms. What can you notice? Are there any significant differences between the “power.0” and “power.1” files?
To save producing thousands of files, a subset can be defined by adding:
>> --bounds $minX $minY $maxX $maxY
Start by viewing a small area, between:
>> minX=507778
>> minY=4267681
>> maxX=509778
>> maxY=4269681
Note that not every file will pass through this area. Those files have been subset to:
/geos/netdata/avtrain/data/3d/active_sensing/week10/gedi_sonoma_area
The script can also output a list of all coordinates in the file.
>> python3 $bin/gediHandler.py --input $filename --writeCoords >> $output
Waveform processing
Before we can use waveform-lidar data, we must process it to remove background noise and identify the ground. We will use the method first proposed in Hofton et al (2000), and described in detail in Wagner et al (2004).
The particular code implementation we are using, “gediMetric”, is described in Hancock et al (2017) and Hancock et al (2015). It is provided on bitbucket, with the README defining the input and output variables.
https://bitbucket.org/StevenHancock/gedisimulator/src/master/
Install with the scripts in the software section. To process a single file, type:
>> gediMetric -input $input -readHDFgedi -ground -varScale 3
Here “varScale” sets the noise threshold based on the background noise statisctics. Type gediMetric –help to see the other denoising options. “–sWidth” and “-minWidth” are particularly useful.
To view the processed waveform, add “-writeFit” on the end. This will generate an ASCII file containing the raw waveform, the denoised waveform and the result of the Gaussian fits. Plot some of these up using either the provided python script (plotFits.py) or your favourite graphing package. The first row defines what each column contains.
Try a few different sets of denoising parameters, being careful to keep track of which outputs are which (keep them in separate folders with different output filename roots) and see how it affects the processed waveform.
The files also contain an estimate of the true ground elevation, from ALS data (note that real GEDI and LVIS data does not). This can be used to assess the ground finding accuracy. A python script is provided to produce a histogram of ground finding accuracy. Use this to find a usable set of denoising parameters (this could be done quickly with a batch script). NOTE when predicting biomass it is more important to minimise the bias than the standard deviation.
>> python3 $bin/groundErrHist.py --input $input --nBins 10
Where “$input” is the metric file (ending .metric.txt) and “$nBins” is the number of bins to give the histogram. What do you notice about the distribution of ground elevation errors? What might be causing them, based on your knowledge of ground finding algorithms.
Once you have a set of denoising parameters you are happy with, a batch script is provided to run over all files.
>> $bin/batchProcess.bash -inDir $inDir/ -outRoot $outRoot -sWidth $sig -varScale $v
Use “-help” to see all options and descriptions. Use “-bounds” to save computational time.
Biomass
Once that is done we can convert to a biomass estimate. Lidar waveform metrics are well correlated (Drake et al 2002). An equation to relate waveform metrics to biomass has been provided by Dr. Duncanson as part of the NASA CMS HiBio project. This was generated by optimising a function to relate ground estimates of biomass to lidar waveform metrics. A scatterplot comparing the ground estimates to the derived lidar estimates gives an idea of the accuracy.
A script is provided for that.
>> python3 $bin/predictBiomass.py --input $input
Where $input is the output of gediMetric. This will produce an ASCII file of “x,y,z,biomass”. It can be displayed in ArcMap, QGIS etc. Note that the biomass estimates here are rough exampled for teaching and should not be taken as truth. The biomass model used is still in development and a final version will be published in the near future.
It can be batch processed to output to a single file:
>> $bin/batchBiomass.bash -inDir $inDir -output $output
Upscaling and further work
We now have footprint level biomass estimates, a vector file. We want use this to generate a raster biomass map. We can take the simple mean.
A script is provided to grid biomass estimates:
>> python3 $bin/gridBiomass.py --input $input --output $output --res $res
Where res is resolution in metres.
This outputs an ASCII raster file. Note that the EPSG is 32610 for displaying in a GIS package.
This assumes that each footprint is an independent sample and that they are representative of what is contained within a given pixel. There are various strategies to remove these assumptions. One approach would be to use a land cover map to stratify the footprint levels estimates, providing a mean biomass value for each land cover type and then weighting those by the relative fraction of that land cover type.
A land cover map could be downloaded from the MODIS data-site, a higher resolution source (eg. https://earthenginepartners.appspot.com/science-2013-global-forest), or we could use the original ALS data to make one ourselves. The raw lidar data over our area of interest is provided in:
This can be run through lasclassify to produce a land cover map. Then the footprint level biomass estimates intersected with it in a GIS package (or python) to get the mean per class per pixel. This is left as an exercise for the student.