1 |
a92ae83f
|
Jim Regetz
|
&args DEMgrid
|
2 |
|
|
|
3 |
|
|
/* Calculates a measure of the noisiness of the surface.
|
4 |
|
|
/* John Gallant 8 Jan 2009
|
5 |
|
|
/*
|
6 |
|
|
/* Modified from noisemag on 15 Oct 2010 to work at 3 second resolution
|
7 |
|
|
/*
|
8 |
|
|
/* The principle is that noise is characterised by short-range variation.
|
9 |
|
|
/* The difference from the local mean will vary erratically over a short
|
10 |
|
|
/* distance in the presence of noise. Terrain features are generally smoother
|
11 |
|
|
/* so do not vary as erratically. The noise is measured as the local standard
|
12 |
|
|
/* deviation of the difference from the mean.
|
13 |
|
|
/*
|
14 |
|
|
/* The local mean is calculated using an annulus so that local peaks and
|
15 |
|
|
/* holes are further accentuated.
|
16 |
|
|
/*
|
17 |
|
|
/* The raw noise magnitude is smoothed using two steps of median filtering
|
18 |
|
|
/* to give the median value over a 25-cell circle.
|
19 |
|
|
/*
|
20 |
|
|
/* The calculated values can be used as the standard deviation of the
|
21 |
|
|
/* noise of the surface.
|
22 |
|
|
/*
|
23 |
|
|
/* One flaw is that higher relief terrain is considered "noisy" by this method.
|
24 |
|
|
/* (Lower relief terrain is not.)
|
25 |
|
|
/*
|
26 |
|
|
/* Could remove the effect of steeper terrain by using a fitted polynomial
|
27 |
|
|
/* rather than a simple mean as the basis for calculating the difference.
|
28 |
|
|
/* Not pursued at this time because the application is to control where
|
29 |
|
|
/* vegetation edges are visible in SRTM data, and high relief areas also
|
30 |
|
|
/* obscure the veg edges, so including the high relief as "high noise"
|
31 |
|
|
/* is not a problem.
|
32 |
|
|
|
33 |
|
|
setwindow %DEMgrid%
|
34 |
|
|
setcell %DEMgrid%
|
35 |
|
|
|
36 |
|
|
noisemean = focalmean(%DEMgrid% , annulus, 1, 2)
|
37 |
|
|
noisediffmn = %DEMgrid% - noisemean
|
38 |
|
|
noiseraw = focalstd(noisediffmn, circle, 2)
|
39 |
|
|
/* setcell minof
|
40 |
|
|
/* noisemag2 = focalmedian(aggregate(noiseraw, 2, median), circle, 2)
|
41 |
|
|
/* setcell %DEMgrid%
|
42 |
|
|
/* noisemag = resample(noisemag2, #, bilinear)
|
43 |
|
|
|
44 |
|
|
/* testing rejection of terrain signal by capturing longer-range differences from mean
|
45 |
|
|
/* noisemag_d3 looks particularly good - terrain effects almost completely absent, for both small and extended features
|
46 |
|
|
|
47 |
|
|
noisemndf = focalmean(noisediffmn, circle, 2) /* calculate mean difference from mean
|
48 |
|
|
noisemndf_abs = sqrt(1 + sqr(noisemndf)) - 1 /* soft absolute value of that
|
49 |
|
|
noisemndf_mdn = focalmedian(noisemndf_abs, circle, 3) /* smooth, this is the terrain signal
|
50 |
|
|
noiseraw_rdc = noiseraw / (1 + 2 * mndiffmn2bfmd) /* reduce noiseraw by the terrain signal
|
51 |
|
|
|
52 |
|
|
setcell minof
|
53 |
|
|
noisemag_2_2 = focalmedian(aggregate(noiseraw_rdc, 2, median), circle, 2)
|
54 |
|
|
setcell %DEMgrid%
|
55 |
|
|
noisemag = resample(noisemag_2_2, #, bilinear)
|