ARL-TR-8271 ● JAN 2018

US Army Research Laboratory

An Automated Energy Detection Algorithm 

Based on Morphological Filter Processing with 

a Semi-Disk Structure

by Kwok F Tom

Approved for public release; distribution is unlimited.
NOTICES

Disclaimers

The findings in this report are not to be construed as an official Department of the 

Army position unless so designated by other authorized documents.

Citation of manufacturer’s or trade names does not constitute an official 

endorsement or approval of the use thereof.

Destroy this report when it is no longer needed. Do not return it to the originator.
ARL-TR-8271 ● JAN 2018

US Army Research Laboratory

An Automated Energy Detection Algorithm 

Based on Morphological Filter Processing with 

a Semi-Disk Structure

by Kwok F Tom 

Sensors and Electron Devices Directorate, ARL

Approved for public release; distribution is unlimited.
ii

REPORT DOCUMENTATION PAGE Form Approved

OMB No. 0704-0188

Public reporting burden for this collection of information is estimated to average 1 hour per response, including the time for reviewing instructions, searching existing data sources, gathering and maintaining the 

data needed, and completing and reviewing the collection information. Send comments regarding this burden estimate or any other aspect of this collection of information, including suggestions for reducing the 

burden, to Department of Defense, Washington Headquarters Services, Directorate for Information Operations and Reports (0704-0188), 1215 Jefferson Davis Highway, Suite 1204, Arlington, VA 22202-4302. 

Respondents should be aware that notwithstanding any other provision of law, no person shall be subject to any penalty for failing to comply with a collection of information if it does not display a currently 

valid OMB control number.

PLEASE DO NOT RETURN YOUR FORM TO THE ABOVE ADDRESS.

1. REPORT DATE (DD-MM-YYYY)

January 2018

2. REPORT TYPE

Technical Report

3. DATES COVERED (From - To)

1 October 2016–30 September 2017

4. TITLE AND SUBTITLE

An Automated Energy Detection Algorithm Based on Morphological Filter 

Processing with a Semi-Disk Structure

5a. CONTRACT NUMBER

5b. GRANT NUMBER

5c. PROGRAM ELEMENT NUMBER

6. AUTHOR(S)

Kwok F Tom

5d. PROJECT NUMBER

5e. TASK NUMBER

5f. WORK UNIT NUMBER

7. PERFORMING ORGANIZATION NAME(S) AND ADDRESS(ES)

US Army Research Laboratory

ATTN: RDRL-SER-E

2800 Powder Mill Road

Adelphi, MD 20783-1138

8. PERFORMING ORGANIZATION REPORT NUMBER

ARL-TR-8271

9. SPONSORING/MONITORING AGENCY NAME(S) AND ADDRESS(ES) 10. SPONSOR/MONITOR'S ACRONYM(S)

11. SPONSOR/MONITOR'S REPORT NUMBER(S)

12. DISTRIBUTION/AVAILABILITY STATEMENT

Approved for public release; distribution is unlimited.

13. SUPPLEMENTARY NOTES

14. ABSTRACT

This report is the result of applying morphological image and statistical processing techniques to the energy detection scenario 

of signals in the RF spectrum domain. A semi-disk structure is applied to the spectral data to estimate the envelope of the 

spectral response in the frequency domain. A morphological “opening” operation is then applied to the spectral envelope to 

obtain some measure of the background noise. These series of techniques provide for an automatic algorithm for determining 

the RF spectral threshold detection level.

15. SUBJECT TERMS

RF spectrum, morphological image processing, detection threshold algorithm, opening, semi-disk structure

16. SECURITY CLASSIFICATION OF:

17. LIMITATION

 OF

 ABSTRACT

UU

18. NUMBER

 OF 

 PAGES

70

19a. NAME OF RESPONSIBLE PERSON

Kwok F Tom

a. REPORT

Unclassified

b. ABSTRACT

Unclassified

c. THIS PAGE

Unclassified

19b. TELEPHONE NUMBER (Include area code)

(301) 394-2612

Standard Form 298 (Rev. 8/98)

Prescribed by ANSI Std. Z39.18
Approved for public release; distribution is unlimited.

iii

Contents

List of Figures v

List of Tables v

Preface vi

1. Introduction 1

2. Data Collection and Statistical Summary 1

3. Statistical Processing 2

3.1 Statistical Analysis 2

3.1.1 Moments 3

3.1.2 Mean 3

3.1.3 Variance 3

3.1.4 Standard Deviation 3

3.1.5 Kurtosis 4

3.1.6 Maximum 4

3.1.7 Minimum 4

3.1.8 Median 4

3.1.9 Rank Order Filter 4

3.1.10 Crest Factor (CF) 5

3.2 Statistical Summary 6

4. Morphological Image Processing 7

5. Algorithm 9

6. Conclusion 11

7. References 12

Appendix A. MATLAB Code 13
Approved for public release; distribution is unlimited.

iv

Appendix B. Graphs of Morphological Processed RF Spectrum Files 25

Appendix C. Graphs of RF Spectrum Files Calculated Detection 

Threshold 43

List of Symbols, Abbreviations, and Acronyms 61

Distribution List 62
Approved for public release; distribution is unlimited.

v

List of Figures

Fig. 1 Dilation of spectral data file with semi-disk structure ........................ 10

Fig. 2 Close-in examination of dilated spectral data file............................... 11

List of Tables

Table 1 Summary of the RF spectrum collection .............................................. 2

Table 2 Statistical analysis of the RF spectrum measurements......................... 7
Approved for public release; distribution is unlimited.

vi

Preface

Energy detection in the RF spectrum is the most basic technique for signal 

detection. Typically, this requires establishing an energy detection threshold based 

on a noise-only condition (i.e., no signal present in the RF spectrum). Initially, the 

area of exploration for this report was to examine the potential of an automatic 

energy detection thresholding algorithm based on the RF measurement. It would 

not require the preliminary RF spectrum noise-only measurements to establish the 

energy detection threshold.

The technique of “opening” was very successful in determining an energy detection 

threshold for the RF spectrum with signal and noise-only environments. After 

examining the results of the energy detection algorithm, a second algorithm that 

incorporated a semi-disk structure on RF spectrum was developed. The windowing 

of the RF spectral data with this semi-disk structure had the benefit of reducing the 

number of false alarms on noise. This is the second of 5 reports that detail the 

energy detection techniques examined with the recorded RF spectrum 

measurements.

1–4

 1 Tom K. An automated energy detection algorithm based on morphological and statistical processing 

techniques. Adelphi (MD): Army Research Laboratory (US); 2018 Jan. Report No.: ARL-TR-8272. 2 Tom K. An automated energy detection algorithm based on morphological filter processing with a

modified watershed transform. Adelphi (MD): Army Research Laboratory (US); 2018 Jan. Report No.: ARLTR-8270. 3 Tom K. An automated energy detection algorithm based on kurtosis-histogram excision. Adelphi (MD): 

Army Research Laboratory (US); 2018 Jan. Report No.: ARL-TR-8269. 4 Tom K. An automated energy detection algorithm based on consecutive mean excision. Adelphi (MD): 

Army Research Laboratory (US); 2018 Jan. Report No.: ARL-TR-8268.
Approved for public release; distribution is unlimited.

1

1. Introduction

Energy detection is the simplest method of detecting a signal in the frequency 

spectrum. When doing so, a comparison is made between the frequency spectral 

component energy and a detection threshold level. Knowledge of the frequency 

spectrum is usually necessary to establish this detection threshold level. Various 

methods can be used to determine the detection threshold level, such as establishing 

the system noise statistics offline and setting the threshold for a given probability 

of detection versus probability of false alarm. 

The goal of this study was to develop an algorithm that can analyze the frequency 

spectrum and establish a threshold value based on the spectral data. The automated 

processing employs techniques from image processing and statistical analysis. A 

combination of techniques from these 2 areas resulted in an algorithm for 

determining a threshold detection level. In this evaluation, a dilation of the spectral 

data was performed with a semi-disk structure. This operation provided an estimate 

of the spectral data envelope in the frequency domain. Morphological filtering was 

then applied to the spectral envelope data to establish a noise level. The 

morphological filtering, consisting of erosion and dilation, is repeatedly applied to 

the spectral envelope data until certain statistical criteria are met. An offset was 

then added to the morphological filtered data to form the threshold for energy 

detection.

2. Data Collection and Statistical Summary

The local RF spectrum was measured in 2013 on the rooftop of building 204 at the 

US Army Research Laboratory’s (ARL’s) Adelphi location. An Agilent N9342CN 

spectrum analyzer and a Discone antenna was used to collect RF spectrum data. 

This spectrum analyzer was operated under the control of a LabVIEW software 

program to acquire and store data with different resolution bandwidth (RBW) from 

1 kHz to 1 MHz.

These data files represent various sizes of RF spectral coverage from 10 MHz up 

to 4 GHz. The number of data files for each spectral band varied. The larger RBW 

files were acquired over seconds of data acquisition time versus the small RBW 

files. Small RBW provides fine spectral resolution, but it impacts data acquisition 

time for spectral coverage and data size. Depending on the RBW and spectral 

coverage, a data file could require a few hours of acquisition.

Table 1 summarizes the data collection measurement files. The RF spectral bands 

covered the spectrum from 10 MHz to 4 GHz. In general, various spectral bands 
Approved for public release; distribution is unlimited.

2

were measured with 4 RBW configuration. Data file size is inversely proportional 

to the RBW. Data file size is proportional to the spectral band coverage. The data 

size varied from approximately 1 KSample to 4 MSample data points per file.

Table 1 Summary of the RF spectrum collection

Spectral band coverage

Number of RBW 

measurement RBW

1.1–1.6 GHz 4 1 kHz, 10 kHz, 100 kHz, 

1 MHz

2–3 GHz 4 1 kHz, 10 kHz, 100 kHz, 

1 MHz

3–4 GHz 3 10 kHz, 100 kHz, 1 MHz

4–6 GHz 3 10 kHz, 100 kHz, 1 MHz

10 MHz–1 GHz 4 1 kHz, 10 kHz, 100 kHz, 

1 MHz

10 MHz–2 GHz 4 1 kHz, 10 kHz, 100 kHz, 

1 MHz

10 MHz–3 GHz 4 1 kHz, 10 kHz, 100 kHz, 

1 MHz

10 MHz–4 GHz 4 1 kHz, 10 kHz, 100 kHz, 

1 MHz

100 MHz–1 GHz 1 100 kHz

3. Statistical Processing

3.1 Statistical Analysis

Statistical analysis is the mathematical science dealing with the analysis or 

interpretation of data. The data analyst uses a few straightforward statistical 

techniques as a means of summarizing the collected data. These statistical 

techniques are under the area of descriptive statistics, which is a methodology to 

condense the data in quantitative terms.

In commercial prognostics and diagnostic vibrational monitoring applications, 

statistical techniques that are mainly used for alarm purposes in industrial plants 

are the statistical moments of order 2, 3, and 4. The probability density function 

(PDF) of the vibrational time series of a good bearing has a Gaussian distribution 

(also known as a normal distribution), whereas a damaged bearing results in a nonGaussian distribution with dominant tails because of a relative increase in the 

number of high levels of acceleration. These techniques can be applied to the RF 

spectral data with a different interpretation of the results.1
Approved for public release; distribution is unlimited.

3

3.1.1 Moments

If these moments are calculated about the mean, they are called central statistical 

moments. The first and second moments are well known, being the mean and the 

variance, respectively. These are analogous to the first and second area moments of 

inertia with the area shape defined by the PDF. The third moment is termed 

skewness and the fourth moment is termed kurtosis. The general equation for the 

order of moment is as follows:

MMpp = 1

NN�(xxii − xx�)pp

NN

ii=1

,

wwheeeeee pp iiii tthee oo oo oooo tthee mm mm mm ,

NN iiii tthee nn nn oooo dd dd vv vv ,

ii iiii tthee ii oooo tthee dd dd vv vv ,

xx� iiii tthee mm mm vv vv oooo tthee dd dd .

3.1.2 Mean

Mean is the most common measure of a statistical distribution. In this case, mean 

is the arithmetic average for a set of measurements. 

xx� = μμ = 1

NN�xxii

NN

ii=1

.

3.1.3 Variance

Variance is a measure of the dispersion of a waveform about its mean—also called 

the second moment of the measurements.

σσ2 = 1

NN�(xxii − xx�)2

NN

ii=1

.

3.1.4 Standard Deviation

Standard deviation is a measure of the variation of a set of data values. The standard 

deviation is defined as the square root of the variance moment.
Approved for public release; distribution is unlimited.

4

σσ = �1

NN�(xxii − xx�)2

NN

ii=1

.

3.1.5 Kurtosis

Kurtosis is the fourth statistical moment, normalized by the standard deviation to 

the fourth power. It is a measure of whether the data are peaked or flat relative to a 

normal distribution. The noise in the RF spectrum is typically considered to have a 

normal distribution. The normal distribution has a value of 3.

κκ = MM4

σσ4 .

κκ =

1

NN ∑ (xxii − xx�) NN 4

ii=1

σσ4 .

κκ = 1

NNσσ4�(xxii − xx�)4

NN

ii=1

.

3.1.6 Maximum

“Max” is the largest value of a set of numbers.

yy = max[xx( )].

3.1.7 Minimum

“Min” is the smallest value of a set of numbers.

yy = min[xx( )].

3.1.8 Median

The statistical median is an order statistic that gives the “middle” value of a set of 

samples. Median is the middle value of a set of data values that divides the set into 

2 groups. Half the groups exist below and half exist above this value.

3.1.9 Rank Order Filter 

Rank order filter is a sorting process by which a set of numbers is ordered from the 

smallest to the largest value. Rank order filtering is a nonlinear filtering technique 

that orders the contents of a filter kernel and selects the sample indexed by rank 

from the magnitude ordered samples.
Approved for public release; distribution is unlimited.

5

Rank order filtering can be summarized as follows:

yy = ∑ iixx�(ii)

NN

ii=1

where xx�(ii), = 1, ⋯ , NN is the result of sorting the data in ascending order. With this 

definition, the max, min, and median can be obtained from a rank order filter as 

follows2:

Min:

yy = ∑ iixx�(ii)

NN

ii=1

ii = 1 ii = 1

0 oo heeee ee

Max:

yy = ∑ iixx�(ii)

NN

ii=1

ii = 1 ii = NN

0 oo heeee ee

Median:

yy = ∑ iixx�(ii)

NN

ii=1

ii = 1 ii = NN

2

0 oo heeee ee

3.1.10 Crest Factor (CF)

Crest factor (CF) is a measure of a waveform showing the ratio of peak values to 

the effective value. In other words, CF indicates how extreme the peaks are in a 

waveform.

CC FF FF = |xx|pppppppp

xxrrrrrr

.

Noise sources are characterized by their CF, which is the peak to average ratio of 

the noise. In a technical bulletin, XiTRON reported CF values between 5 and 7 for 

random noise.3 For example, a 5:1 CF of the noise voltage is 20log5 = 14 dB. This 

is a measure of the quality of the noise distributions and one way to measure its 
Approved for public release; distribution is unlimited.

6

Gaussian nature.4 For the purpose of algorithm development, the CF equation was 

modified as follows:

CC FF FF = MM

MM .

3.2 Statistical Summary 

There were 31 different groupings of the RF spectrum data measurements. The 

number of data files under each of the main groupings was varied. For the purpose 

of developing the algorithm, only a single data file was selected from each group. 

Each data file was processed to obtain the following characteristics: RBW, mean, 

standard deviation, median, max, min, kurtosis, and CF. A summary of the results 

is in Table 2.

The following results were noted:

• The smaller the RBW, the lower the noise floor. The equation for thermal 

noise power is PP = kkkkkk, where B is bandwidth. In this case, the RBW of 1 

kHz has the lowest noise value. The RBW of 1 MHz has the highest noise 

value.

• Each 10-fold increase in bandwidth results in a 10-dB increase in noise 

power. This relationship is illustrated in this data set.

• The calculated mean and median values are very close for a given RF 

measurement configuration.

• The CF for noise-only data files was on the order of approximately 10 to 13. 

Noise-only data files were estimated by visually inspecting the spectrum plot.
Approved for public release; distribution is unlimited.

7

Table 1 Statistical analysis of the RF spectrum measurements

4. Morphological Image Processing

A description of a technique for automatically estimating the noise floor spectrum 

was given in a conference back in 1997.5 This technique is based on applying the 

morphological binary image processing operators to the RF spectrum. The 

technique works well in both flat and nonflat noise floor spectra. Morphology 

image processing is a set of nonlinear operations. Ready et al. note that “humans 

are good at estimating the noise floor spectrum by ‘eyeballing’ a spectral plot. 

Intuitively, we separate the spectral humps from the noise floor spectrum by 

eliminating those parts of the spectrum shape that are due to signals and visually 

draw in the noise floor spectrum.”5

Morphology is a broad set of processing techniques that process images based on 

shapes. Morphological operations apply a structuring element to an input image, 

creating an output image of the same size. In a morphological operation, the value 

of each pixel in the output image is based on a comparison of the corresponding 

pixel in the input image with its neighbors. By choosing the size and shape of the 

Filename Spectral Band RBW Mean SD Median Max Min Kurtosis CF

Air_test_1.1GHz_1.6GHzb_03_28_14_06_33_22 1.1 - 1.6 GHz 1000 -112.694 3.8444 -112.3 -101.2 -138.4 3.5571 11.1

Air_test_1.1GHz_1.6GHza_03_27_14_07_14_05 1.1 - 1.6 GHz 10000 -102.12 3.7226 -101.8 -91.69 -124.5 3.5134 10.11

Air_test_1.1GHz_1.6GHzc_03_31_14_06_47_26 1.1 - 1.6 GHz 100000 -90.9109 2.597 -90.75 -83.13 -101.2 3.0217 7.62

Air_test_1.1GHz_1.6GHzd_04_01_14_06_54_03 1.1 - 1.6 GHz 1000000 -80.5634 3.5741 -80.745 -55.36 -89.8 14.6861 25.385

Air_test_2GHz_3GHza_06_05_14_07_15_58 2 - 3 GHz 1000 -111.219 4.4991 -111.1 -77.55 -139 6.8918 33.55

Air_test_2GHz_3GHzb_05_29_14_06_28_28 2 - 3 GHz 10000 -101.006 4.1808 -100.8 -74.67 -122 6.1995 26.13

Air_test_2GHz_3GHzc_05_29_14_04_09_27 2 - 3 GHz 100000 -89.5666 3.4154 -89.68 -64.27 -100.3 13.2282 25.41

Air_test_2GHz_3GHzd_05_28_14_00_06_18 2 - 3 GHz 1000000 -79.6014 3.0753 -79.6 -58.74 -88.88 8.4466 20.86

Air_test_3GHz_4GHzb_06_12_14_06_30_17 3 - 4 GHz 10000 -100.51 3.7072 -100.2 -89.68 -123.7 3.5901 10.52

Air_test_3GHz_4GHzC_06_11_14_06_44_40 3 - 4 GHz 100000 -89.1625 2.5513 -89.075 -80.68 -102.3 3.0774 8.395

Air_test_3GHz_4GHzd_06_09_14_07_09_51 3 - 4 GHz 1000000 -79.0828 2.445 -78.94 -72.76 -88.5 2.9651 6.18

Air_test_4GHz_6GHza_04_24_14_07_05_26 4 - 6 GHz 1000 -107.559 4.1793 -107.2 -94.31 -136.3 3.3696 12.89

Air_test_4GHz_6GHzb_04_28_14_06_23_12 4 - 6 GHz 10000 -97.245 4.0398 -96.95 -84.05 -119.7 3.3249 12.9

Air_test_4GHz_6GHzc_04_29_14_06_48_42 4 - 6 GHz 100000 -85.9125 3.082 -85.77 -75.84 -98.78 2.9068 9.93

Air_test_10MHz_1GHza_04_17_14_07_09_30 10 MHz - 1 GHz 1000 -110.348 8.6126 -112.1 -28.37 -140.3 6.4788 83.73

Air_test_10MHz_1GHzb_04_21_14_07_17_44 10 MHz - 1 GHz 10000 -100.168 8.4775 -102 -27.95 -126.7 6.7602 74.05

Air_test_10MHz_1GHzc_04_22_14_06_39_58 10 MHz - 1 GHz 100000 -87.948 8.8113 -90.74 -27.28 -103.6 6.513 63.46

Air_test_10MHz_1GHzd_04_23_14_06_36_16 10 MHz - 1 GHz 1000000 -78.1054 8.3327 -80.59 -27.82 -89.88 7.5282 52.77

Air_test_10MHz_2GHzd_03_06_14_08_52_04 10 MHz - 2 GHz 1000 -111.275 6.9798 -112.1 -27.34 -138.3 9.5085 84.76

Air_test_10MHz_2GHzC_02_27_14_06_43_55 10 MHz - 2 GHz 10000 -100.894 7.0157 -101.8 -26.58 -125.9 10.4224 75.22

Air_test_10MHz_2GHzB_02_25_14_06_35_53 10 MHz - 2 GHz 100000 -89.242 6.7487 -90.59 -26.09 -102.4 12.7961 64.5

Air_test_10MHz_2GHz_02_20_14_06_55_36 10 MHz - 2 GHz 1000000 -78.8448 7.3203 -80.6 -26.54 -92.41 11.3819 54.06

Air_test_10MHz_3GHz_03_11_14_06_32_51 10 MHz - 3 GHz 1000 -111.373 6.1356 -111.8 -29.42 -140.5 11.0474 82.38

Air_test_10MHz_3GHzB_03_13_14_06_32_45 10 MHz - 3 GHz 10000 -100.787 6.3591 -101.4 -25.36 -125.9 11.7254 76.04

Air_test_10MHz_3GHzC_03_18_14_06_32_08 10 MHz - 3 GHz 100000 -89.1614 6.0257 -90.14 -24.16 -102.4 15.3376 65.98

Air_test_10MHz_3GHzd_03_20_14_08_09_07 10 Mhz - 3 GHz 1000000 -78.6492 6.7418 -80.8 -25.81 -90.38 13.0673 54.99

Air_test_10MHz_4GHza_04_10_14_07_01_22 10 MHz - 4 GHz 1000 -111.158 5.6813 -111.4 -28.22 -143.1 10.9573 83.18

Air_test_10MHz_4GHzb_04_11_14_06_10_18 10 MHz - 4 GHz 10000 -100.838 5.5851 -101.1 -26.85 -128 12.2576 74.25

Air_test_10MHz_4GHzc_04_15_14_07_01_09 10 MHz - 4 GHz 100000 -89.4431 5.1076 -90.03 -25.05 -103.5 18.317 64.98

Air_test_10MHz_4GHzd_04_16_14_06_42_27 10 MHz - 4 GHz 1000000 -79.0396 5.5441 -79.89 -26.72 -91.64 16.9073 53.17

Air_test_100MHz_1GHz_02_19_14_07_20_46 100 MHz - 1 GHz 100000 -88.0772 9.0701 -90.82 -23.91 -102.2 6.994 66.91
Approved for public release; distribution is unlimited.

8

neighborhood, you can construct a morphological operation that is sensitive to 

specific shapes in the input image.

The most basic morphological operations are dilation and erosion. Dilation adds 

pixels to the boundaries of the objects in an image, while erosion removes pixels 

on object boundaries. The number of pixels added or removed from the objects in 

an image depends on the size and shape of the structuring element used to process 

the image. In the morphological dilation and erosion operations, the state of any 

given pixel in the output image is determined by applying a rule to the 

corresponding pixel and its neighbors in the input image. These rules are known as 

the erosion and dilation:

Erosion:

AA ⊖ = ⋂

∈

AA−bb .

Dilation:

AA⨁ = ⋃

∈ AA .

The opening of A by B is obtained by the dilation of A by B, followed by erosion of 

the resulting structure by B.

Opening:

AA ∘ = (AA ⊖ )⨁ .

The application of the opening technique on the RF spectrum is the basis for noise 

floor estimation. The estimation of the noise floor is used as a reference level to set 

a threshold detection level.6

In signal, statistical, and image processing, the minimum and maximum operators 

are typically encountered. Gil and Kimmel note, “In mathematical morphology, the 

result of such an operator is referred to as the erosion (or dilation) of the signal with 

a structuring element given by a pulse of width p.”7

 

For the 1-D case, this reduces to a simple filter of just providing the max or min 

value of a set of values.7

1-D max filter: Given a sequence xx0, ⋯ , xxnn−1 and an integer pp > 1, compute 

yyii = max

0≤jj<pp

xxii+jj

for = 0, ⋯ , − pp .
Approved for public release; distribution is unlimited.

9

1-D min filter: Given a sequence xx0, ⋯ , xxnn−1 and an integer pp > 1, compute

yyii = min

0≤jj<pp

xxii+jj

for = 0, ⋯ , − pp .

5. Algorithm

The algorithm development process was executed in MATLAB. Appendix A is the 

code used to process and generate the enclosed RF spectrum signature with the 

corresponding results of the morphological processing and resultant threshold 

detection level. The following is a description of the thresholding detection level 

generation based on morphological filter processing:

1) Determine the RBW of the spectral data.

2) Determine some statistics on the spectral data file: median, max, and CF.

3) Form an array based on a semi-disk. This array represents an 8-segment 

structure of a semi-disk at an arc of 22.5°.

4) Using the semi-disk structure, dilate the spectral data file. At each point in 

the origin spectral data, the data array of the semi-disk is added to the 

spectral data, and the maximum value saved. Figure 1 is a graph of a spectral 

data file with the dilated data overlaid on the spectral data file. The blue 

color is the measured spectral data, and the red color is the dilated results. 

Note how the dilated file follows the maximum envelop of the spectral data. 

Figure 2 is a zoomed-in section of this data file. One can see that valleys 

tend to be filled, and the RF spectral data take on a repeated semi-disk 

pattern.

5) Perform an erosion operation on the semi-disk dilated spectral data as 

computed in step 4. Starting from the lowest-frequency component, create 

the eroded spectral array with an initial sliding window of size 2. Basically, 

calculate the minimum value for the given window size, fill in that position,

and shift the processing window one frequency position to calculate the next 

value. Continue on across the entire semi-disk spectral data array. For the 

last value, use the previous valid value that is calculated.

6) Perform a dilation operation on the eroded spectral array as calculated in 

step 5. In this case, the calculations are performed starting at the highest 

frequency and continue to the lowest frequency. Dilation is simply the 

maximum value in the windowed data set. The initial window size is 2. In 
Approved for public release; distribution is unlimited.

10

this case, the maximum value is determined from the sliding window data 

and used to create the dilated spectral array.

7) Calculate statistics on the dilated spectral array from step 6: mean, median, 

max, and CF.

8) Repeat the morphological filtering of erosion and dilation on the spectral 

array by increasing the window size by 1 for every pass. For example, on 

the second pass, the window size is 3 for erosion and dilation operations. 

Continue with the morphological filtering process as long as the CF is 

greater than 10 and the mean power greater than 0.01 on each pass of the 

morphological filtering.

9) If the CF is not greater than 10, then perform the morphological filtering

operations until the average power is greater than 0.01.

10) Once the convergence has been met, use the resulting morphological 

filtered data array (i.e., the dilated spectral array) to form the threshold as 

follows: 

h hoo oo =

MM hoo oooo FF SS AA + 25 −

2.9 x 10 log10 �

RRRRRR

1000�

Fig. 1 Dilation of spectral data file with semi-disk structure
Approved for public release; distribution is unlimited.

11

Fig. 2 Close-in examination of dilated spectral data file

6. Conclusion

It is possible to estimate the energy detection threshold level based on 

morphological filtering. 

Appendix B displays the graphs of the morphological processing for each of the 

selected RF spectrum data files. The red curve is the result of the morphological 

processing overlaid on top of the RF spectrum signature. 

Appendix C displays the results of the threshold overlaid on top of the RF spectrum 

signature. The graphs are intended to provide the reader a qualitative sense of the 

effectiveness of the automatic enhanced threshold generation algorithm. Reduction 

in false alarms (over just morphological filtering case) can be achieved with an 

initial dilation of the spectral data with a semi-disk structure. Actual processing 

time and the number of iterations through the morphological filtering routine are 

highly dependent on the RBW and spectral coverage band. A small RBW over a 

large coverage band requires more time to complete the process. The threshold 

detection level obtained works well for signals present or noise-only RF spectrum 

data files. The red curve is the threshold when added to the morphological 

processed array as shown in Appendix B. 

Overall, a visual inspection of the graphs shows that the algorithm works well. 

Further refinement of the algorithm is necessary to address potential issues related 

to the boundary edge of the spectral data file. The stop criteria for the morphological 

filtering and thresholding equation can be optimized further to reduce the false 

alarm detections.
Approved for public release; distribution is unlimited.

12

7. References

1. Tom KF. A primer on vibrational ball bearing feature generation for 

prognostics and diagnostics algorithms. Aberdeen Proving Ground (MD): 

Army Research Laboratory (US); 2015 Mar. Report No.: ARL-TR-7230.

2. Mitra AK, Lewis, TL, Paul AS. Arnab K. Ultra-wideband radar data models 

and target detection with adaptive rank-order filters. Proc SPIE 4727; 2002.

3. Crest factors and power analyzers. San Diego (CA): XiTRON Technologies; 

nd. Technical Brief TB 103 [accessed 2017 Nov 8]. 

http://www.xitrontech.com/assets/002/5788.pdf.

4. Noise: frequently asked questions. East Hanover (NJ): Noisewave; nd 

[accessed 2017 Nov 8]. http://www.noisewave.com/faq.pdf.

5. Ready MJ, Downey ML, Corbalis LJ. Automatic noise floor spectrum 

estimation in the presence of signals. Proceedings of Thirty-First Asilomar 

Conference on Signals, Systems & Computers; 1997 Nov. p 877–881.

6. Sequeira S. Energy based spectrum sensing for enabling dynamic spectrum 

access in cognitive radios [master’s thesis]. [Camden (NJ)]: Rutgers 

University; 2011 May.

7. Gil Y, Kimmel R. Efficient dilation, erosion, opening, and closing algorithms.

IEEE Transactions on Pattern Analysis and Machine Intelligence.

2002;24(12): 1606–1617.
Approved for public release; distribution is unlimited.

13

Appendix A. MATLAB Code
Approved for public release; distribution is unlimited.

14

function DiskMorphologicalFilterDetection5report()

% Make selection for data file to process

% SELECT = 1 to 31

SELECT = 1;

dir{1} = 'K:\CognitiveRadar\spectrum monitoring\data\building204-

4c085\Air_test_1.1GHz_1.6GHza\';

filename{1} = 'Air_test_1.1GHz_1.6GHza_03_27_14_07_14_05';

dir{2} = 'K:\CognitiveRadar\spectrum monitoring\data\building204-

4c085\Air_test_1.1GHz_1.6GHzb\';

filename{2} = 'Air_test_1.1GHz_1.6GHzb_03_28_14_06_33_22';

dir{3} = 'K:\CognitiveRadar\spectrum monitoring\data\building204-

4c085\Air_test_1.1GHz_1.6GHzc\';

filename{3} = 'Air_test_1.1GHz_1.6GHzc_03_31_14_06_47_26';

dir{4} = 'K:\CognitiveRadar\spectrum monitoring\data\building204-

4c085\Air_test_1.1GHz_1.6GHzd\';

filename{4} = 'Air_test_1.1GHz_1.6GHzd_04_01_14_06_54_03';

dir{5} = 'K:\CognitiveRadar\spectrum monitoring\data\building204-

4c085\Air_test_2GHz_3GHza\';

filename{5} = 'Air_test_2GHz_3GHza_06_05_14_07_15_58';

dir{6} = 'K:\CognitiveRadar\spectrum monitoring\data\building204-

4c085\Air_test_2GHz_3GHzb\';

filename{6} = 'Air_test_2GHz_3GHzb_05_29_14_06_28_28';

dir{7} = 'K:\CognitiveRadar\spectrum monitoring\data\building204-

4c085\Air_test_2GHz_3GHzc\';

filename{7} = 'Air_test_2GHz_3GHzc_05_29_14_04_09_27';

dir{8} = 'K:\CognitiveRadar\spectrum monitoring\data\building204-

4c085\Air_test_2GHz_3GHzd\';

filename{8} = 'Air_test_2GHz_3GHzd_05_28_14_00_06_18';

dir{9} = 'K:\CognitiveRadar\spectrum monitoring\data\building204-

4c085\Air_test_3GHz_4GHzb\';

filename{9} = 'Air_test_3GHz_4GHzb_06_12_14_06_30_17';

dir{10} = 'K:\CognitiveRadar\spectrum 

monitoring\data\building204-4c085\Air_test_3GHz_4GHzC\';

filename{10} = 'Air_test_3GHz_4GHzC_06_11_14_06_44_40';

dir{11} = 'K:\CognitiveRadar\spectrum 

monitoring\data\building204-4c085\Air_test_3GHz_4GHzd\';

filename{11} = 'Air_test_3GHz_4GHzd_06_09_14_07_09_51';
Approved for public release; distribution is unlimited.

15

dir{12} = 'K:\CognitiveRadar\spectrum 

monitoring\data\building204-4c085\Air_test_4GHz_6GHza\';

filename{12} = 'Air_test_4GHz_6GHza_04_24_14_07_05_26';

dir{13} = 'K:\CognitiveRadar\spectrum 

monitoring\data\building204-4c085\Air_test_4GHz_6GHzb\';

filename{13} = 'Air_test_4GHz_6GHzb_04_28_14_06_23_12';

dir{14} = 'K:\CognitiveRadar\spectrum 

monitoring\data\building204-4c085\Air_test_4GHz_6GHzc\';

filename{14} = 'Air_test_4GHz_6GHzc_04_29_14_06_48_42';

dir{15} = 'K:\CognitiveRadar\spectrum 

monitoring\data\building204-4c085\Air_test_10MHz_1GHza\';

filename{15} = 'Air_test_10MHz_1GHza_04_17_14_07_09_30';

dir{16} = 'K:\CognitiveRadar\spectrum 

monitoring\data\building204-4c085\Air_test_10MHz_1GHzb\';

filename{16} = 'Air_test_10MHz_1GHzb_04_21_14_07_17_44';

dir{17} = 'K:\CognitiveRadar\spectrum 

monitoring\data\building204-4c085\Air_test_10MHz_1GHzc\';

filename{17} = 'Air_test_10MHz_1GHzc_04_22_14_06_39_58';

dir{18} = 'K:\CognitiveRadar\spectrum 

monitoring\data\building204-4c085\Air_test_10MHz_1GHzd\';

filename{18} = 'Air_test_10MHz_1GHzd_04_23_14_06_36_16';

dir{19} = 'K:\CognitiveRadar\spectrum 

monitoring\data\building204-4c085\Air_test_10MHz_2GHz\';

filename{19} = 'Air_test_10MHz_2GHz_02_20_14_06_55_36';

dir{20} = 'K:\CognitiveRadar\spectrum 

monitoring\data\building204-4c085\Air_test_10MHz_2GHzB\';

filename{20} = 'Air_test_10MHz_2GHzB_02_25_14_06_35_53';

dir{21} = 'K:\CognitiveRadar\spectrum 

monitoring\data\building204-4c085\Air_test_10MHz_2GHzC\';

filename{21} = 'Air_test_10MHz_2GHzC_02_27_14_06_43_55';

dir{22} = 'K:\CognitiveRadar\spectrum 

monitoring\data\building204-4c085\Air_test_10MHz_2GHzd\';

filename{22} = 'Air_test_10MHz_2GHzd_03_06_14_08_52_04';

dir{23} = 'K:\CognitiveRadar\spectrum 

monitoring\data\building204-4c085\Air_test_10MHz_3GHz\';

filename{23} = 'Air_test_10MHz_3GHz_03_11_14_06_32_51';

dir{24} = 'K:\CognitiveRadar\spectrum 

monitoring\data\building204-4c085\Air_test_10MHz_3GHzB\';

filename{24} = 'Air_test_10MHz_3GHzB_03_13_14_06_32_45';
Approved for public release; distribution is unlimited.

16

dir{25} = 'K:\CognitiveRadar\spectrum 

monitoring\data\building204-4c085\Air_test_10MHz_3GHzC\';

filename{25} = 'Air_test_10MHz_3GHzC_03_18_14_06_32_08';

dir{26} = 'K:\CognitiveRadar\spectrum 

monitoring\data\building204-4c085\Air_test_10MHz_3GHzd\';

filename{26} = 'Air_test_10MHz_3GHzd_03_20_14_08_09_07';

dir{27} = 'K:\CognitiveRadar\spectrum 

monitoring\data\building204-4c085\Air_test_10MHz_4GHza\';

filename{27} = 'Air_test_10MHz_4GHza_04_10_14_07_01_22';

dir{28} = 'K:\CognitiveRadar\spectrum 

monitoring\data\building204-4c085\Air_test_10MHz_4GHzb\';

filename{28} = 'Air_test_10MHz_4GHzb_04_11_14_06_10_18';

dir{29} = 'K:\CognitiveRadar\spectrum 

monitoring\data\building204-4c085\Air_test_10MHz_4GHzc\';

filename{29} = 'Air_test_10MHz_4GHzc_04_15_14_07_01_09';

dir{30} = 'K:\CognitiveRadar\spectrum 

monitoring\data\building204-4c085\Air_test_10MHz_4GHzd\';

filename{30} = 'Air_test_10MHz_4GHzd_04_16_14_06_42_27';

dir{31} = 'K:\CognitiveRadar\spectrum 

monitoring\data\building204-4c085\Air_test_100MHz_1GHz\';

filename{31} = 'Air_test_100MHz_1GHz_02_19_14_07_20_46';

str_meta = 

sprintf('%s%s.mspecrawdata',dir{SELECT},filename{SELECT});

str_data = 

sprintf('%s%s.specrawdata',dir{SELECT},filename{SELECT});

% Get metadata on selected data file

fid_meta = fopen(str_meta);

META = textscan(fid_meta,'%s');

ave = META{1}{1};

ref = META{1}{2};

attn = META{1}{3};

rbw = META{1}{4};

% Read in data file

fid_data = fopen(str_data);

g=0;

f=0;

a=0;

while(g==0)

 ft=fgetl(fid_data);

 f=[f,str2num(ft)];

 at=fgetl(fid_data);
Approved for public release; distribution is unlimited.

17

 a=[a,str2num(at)];

 g=feof(fid_data);

end

f=f(2:end);

a=a(2:end);

InputIndex = length(a);

display (InputIndex);

InputArray = a;

% Determine resolution bandwidth

Bandwidth = f(3) - f(2);

% Plot data file

figure(1);

plot(f/1e6,a);

axis tight;

grid;

xlabel('Frequency (MHz)');

ylabel('Power (dBm)');

title(strcat(rbw,'--',ave,'--',ref,'--

',attn),'Interpreter','none');

% Determine some statistics values for data

M = mean(a);

Med = median(a);

S = std(a);

Max = max(a);

Min = min(a);

Kurt = kurtosis(a);

Range = abs(Min - Max);

Number = floor(Range);

Bins = Number * 2;

disp([M S Med Max Min Kurt]);

% Generate semi-disk structure

 part = 8;

 radius = 3;

 

 for i = 0: part
Approved for public release; distribution is unlimited.

18

 ang = i*22.5; 

 x = radius * cosd(ang);

 y = radius * sind(ang);

 xarray(i+1) = x;

 yarray(i+1) = y;

 end

 

 diskarray = -1*yarray;

for m = 1:InputIndex

 InputArray(m+part/2) = a(m); 

end

for m = 1:part/2

 InputArray(m) = Min; 

end

 

NewInputArray = length(InputArray);

for m = 1:part/2

 InputArray(NewInputArray+m) = Min;

end

 

 

TempTopHatArray = InputArray;

j=part;

% Apply semi-disk structure to RF spectrum

 for k = 1:InputIndex

 arraysum = InputArray(k:k+j)+ diskarray;

 TempTopHatArray(k) = max(arraysum);

 end

NewInputArray = length(TempTopHatArray);

 

display(NewInputArray)

for m = 1:InputIndex

 TopHatArray(m) = TempTopHatArray(m+part/2); 

end
Approved for public release; distribution is unlimited.

19

CountTopHat = length(TopHatArray);

display(CountTopHat);

 

% Plot RF spectrum with semi-disk structure applied

 figure(20);

 plot(TopHatArray);

% Plot orignal RF spectrum with overlay of the RF spectrum with 

the

% semi-disk structure applied

figure(30);

plot(f/1e6,a);

axis tight;

grid;

xlabel('Frequency (MHz)');

ylabel('Power (dBm)');

title(strcat(rbw,'--',ave,'--',ref,'--

',attn),'Interpreter','none');

hold on

plot(f/1e6,TopHatArray,'r');

hold off;

 

ErosionIndex = 1;

ErosionArray = a;

DilationArray = a;

InputArray = TopHatArray;

% Plot RF spectrum with semi-disk structure applied

figure(2);
Approved for public release; distribution is unlimited.

20

plot(InputArray);

PreviousPower = sum(a)/InputIndex;

display (PreviousPower);

% Initialize parameters

j= 1;

DiffPower = 1;

DilationMedian = Med;

DilationMax = Max;

DilationCF = DilationMax-DilationMedian;

 

if DilationCF > 10 && j == 1

 display('CF'); 

% Perform Erosion operation on data file

 while DilationCF > 10

 

 for k = 1:InputIndex-j

 ErosionArray(k) = min(InputArray(k:k+j));

 end

 for k = InputIndex-j:InputIndex

 ErosionArray(k) = InputArray(InputIndex-j);

 end

 

 

% Plot results of Erosion operation 

 figure(5);

 plot(ErosionArray);

% Perform Dilation operation on data file

 DilationArray = ErosionArray;

 for k = InputIndex:j+1

 DilationArray(k) = max(ErosionArray(k-j:k));
Approved for public release; distribution is unlimited.

21

 end

 

 for k = j:1

 DilationArray(k) = ErosionArray(j+1);

 end

 

% Plot results of Dilation operation

 InputArray = DilationArray;

 figure(6);

 plot(DilationArray); 

 

 DilationPower = sum(DilationArray)/InputIndex;

 DilationMedian = median(DilationArray);

 DilationMax = max(DilationArray);

 DiffPower = abs(PreviousPower - DilationPower);

 DilationCF = DilationMax-DilationMedian;

 display ([j DilationPower DilationMedian DilationMax 

DilationCF DiffPower]);

 PreviousPower = DilationPower;

 j= j+1;

 end

 % Check performance of the Erosion & Dilation to see if 

operations 

 % resulted have converged to accept criteria (If power 

difference in 

 % the RF spectrum is less than 0.1 between iterations previous 

and 

 % current processing of RF spectrum

 

 while DiffPower > 0.1 

 

 display('CF & Power');

 

 for k = 1:InputIndex-j

 ErosionArray(k) = min(InputArray(k:k+j));

 end

 for k = InputIndex-j:InputIndex

 ErosionArray(k) = min(InputArray(k:InputIndex));

 end
Approved for public release; distribution is unlimited.

22

 figure(5);

 plot(ErosionArray);

 DilationArray = ErosionArray;

 for k = InputIndex:j+1

 DilationArray(k) = max(ErosionArray(k-j:k));

 end

 

 for k = j:1

 DilationArray(k) = max(ErosionArray(k:1));

 end

 

 InputArray = DilationArray;

 figure(6);

 plot(DilationArray); 

 

 DilationPower = sum(DilationArray)/InputIndex;

 DilationMedian = median(DilationArray);

 DilationMax = max(DilationArray);

 DiffPower = abs(PreviousPower - DilationPower);

 DilationCF = DilationMax-DilationMedian;

 display ([j DilationPower DilationMedian DilationMax 

DilationCF DiffPower]);

 PreviousPower = DilationPower;

 j= j+1;

 

 end

else

 while DiffPower > 0.1 

 display('Power');

 

 for k = 1:InputIndex-j

 ErosionArray(k) = min(InputArray(k:k+j));

 end

 for k = InputIndex-j:InputIndex

 ErosionArray(k) = min(InputArray(k:InputIndex));
Approved for public release; distribution is unlimited.

23

 end

 figure(5);

 plot(ErosionArray);

 DilationArray = ErosionArray;

 for k = InputIndex:j+1

 DilationArray(k) = max(ErosionArray(k-j:k));

 end

 

 for k = j:1

 DilationArray(k) = max(ErosionArray(k:1));

 end

 

 InputArray = DilationArray;

 figure(6);

 plot(DilationArray); 

 

 DilationPower = sum(DilationArray)/InputIndex;

 DilationMedian = median(DilationArray);

 DilationMax = max(DilationArray);

 DiffPower = abs(PreviousPower - DilationPower);

 DilationCF = DilationMax-DilationMedian;

 display ([j DilationPower DilationMedian DilationMax 

DilationCF DiffPower]);

 PreviousPower = DilationPower;

 j= j+1;

 

 end

 

end

% Plot original RF spectrum and threshold

figure(3);

plot(f/1e6,a);

axis tight;

grid;

xlabel('Frequency (MHz)');

ylabel('Power (dBm)');
Approved for public release; distribution is unlimited.

24

title(strcat(filename(SELECT),' -- ',rbw,' 

Hz'),'Interpreter','none');

hold on

% Determine threshold array

Threshold = DilationArray+25-(2.9*log10(Bandwidth));

plot(f/1e6,Threshold,'r');

hold off;

% Plot original RF spectrum and final results of Erosion & 

Dilation

% processing

figure(10);

plot(f/1e6,a);

axis tight;

grid;

xlabel('Frequency (MHz)');

ylabel('Power (dBm)');

title(strcat(filename(SELECT),' -- ',rbw,' 

Hz'),'Interpreter','none');

hold on

plot(f/1e6,DilationArray,'y');

hold off;

% Save results to data file

lab = num2str(SELECT);

label = strcat('Disk2-Figure',lab,'.jpg');

saveas(gcf,label);

end
Approved for public release; distribution is unlimited.

25

Appendix B. Graphs of Morphological Processed RF 

Spectrum Files
Approved for public release; distribution is unlimited.

26
Approved for public release; distribution is unlimited.

27
Approved for public release; distribution is unlimited.

28
Approved for public release; distribution is unlimited.

29
Approved for public release; distribution is unlimited.

30
Approved for public release; distribution is unlimited.

31
Approved for public release; distribution is unlimited.

32
Approved for public release; distribution is unlimited.

33
Approved for public release; distribution is unlimited.

34
Approved for public release; distribution is unlimited.

35
Approved for public release; distribution is unlimited.

36
Approved for public release; distribution is unlimited.

37
Approved for public release; distribution is unlimited.

38
Approved for public release; distribution is unlimited.

39
Approved for public release; distribution is unlimited.

40
Approved for public release; distribution is unlimited.

41
Approved for public release; distribution is unlimited.

42

INTENTIONALLY LEFT BLANK.
Approved for public release; distribution is unlimited.

43

Appendix C. Graphs of RF Spectrum Files Calculated Detection 

Threshold
Approved for public release; distribution is unlimited.

44
Approved for public release; distribution is unlimited.

45
Approved for public release; distribution is unlimited.

46
Approved for public release; distribution is unlimited.

47
Approved for public release; distribution is unlimited.

48
Approved for public release; distribution is unlimited.

49
Approved for public release; distribution is unlimited.

50
Approved for public release; distribution is unlimited.

51
Approved for public release; distribution is unlimited.

52
Approved for public release; distribution is unlimited.

53
Approved for public release; distribution is unlimited.

54
Approved for public release; distribution is unlimited.

55
Approved for public release; distribution is unlimited.

56
Approved for public release; distribution is unlimited.

57
Approved for public release; distribution is unlimited.

58
Approved for public release; distribution is unlimited.

59
Approved for public release; distribution is unlimited.

60
Approved for public release; distribution is unlimited.

61

List of Symbols, Abbreviations, and Acronyms

ARL US Army Research Laboratory

CF crest factor

GHz gigahertz

kHz kilohertz

MHz megahertz

PDF probability density function

RBW resolution bandwidth

RF radio frequency
Approved for public release; distribution is unlimited.

62

1 DEFENSE TECHNICAL

(PDF) INFORMATION CTR

DTIC OCA

2 DIR ARL

(PDF) IMAL HRA 

RECORDS MGMT

RDRL DCL

TECH LIB

1 GOVT PRINTG OFC

(PDF) A MALHOTRA

7 ARL

(PDF) RDRL SER E

M CONN

R DEL ROSARIO

K F TOM

D WASHINGTON

RDRL SER M

J SILVIOUS

RDRL SER U

A MARTONE

RDRL SER W

K RANNEY