10226 FACTA UNIVERSITATIS Series: Electronics and Energetics Vol. 35, No 3, September 2022, pp. 379-391 https://doi.org/10.2298/FUEE2203379L © 2022 by University of Niš, Serbia | Creative Commons License: CC BY-NC-ND Original scientific paper VERIFICATION OF CALCULATION METHOD FOR DRONE MICRO-DOPPLER SIGNATURE ESTIMATION Aleksandar Lebl, Mladen Mileusnić, Dragan Mitić, Jovan Radivojević, Vladimir Matić IRITEL a.d., Belgrade, Serbia Abstract. Drones micro-Doppler signatures obtained by FMCW radars are an excellent procedure for malicious drone detection, identification and classification. There are a number of contributions dealing with recorded spectrograms with these micro-Doppler signatures, but very low number of them has analyzed possibility to calculate echo caused by drone moving parts. In this paper, starting from already existing mathematical apparatus, we presented such spectrograms as a function of changing drone moving parts characteristics: rotor number, blades number, blade length and rotor moving speed. This development is the part of a wider project intended to prevent malicious drone usage. Key words: Malicious drone detection, FMCW radar, Spectrogram, Drone micro- Doppler signatures, calculation method 1. INTRODUCTION Drones or unmanned aerial vehicles (UAVs) fulfil our everyday lives more and more. They may be used in many friendly types of missions as, for example, aerial photography, traffic supervision, disaster monitoring, precise agriculture, industrial inspection, goods delivery and so on. But, on the other side, drones are used for a number of different malicious purposes [1]. Drones may carry explosive devices with the intention to cause numerous victims and damages on objects such as airports, stadiums, governmental buildings, residential areas, commercial and industrial facilities, power plants, etc. They may be used for smuggling activities over state borders or into and out of the prisons, for causing fire in hardly accessible forest areas or to perform assassination on the important persons. There are a huge number of examples for each of these malicious activities, successfully or unsuccessfully realized. This is the reason why systems for drones detection, identification, localization and classification (DILC) become very important today. Received November 25, 2021; revised December 28, 2021; accepted January 10, 2022 Corresponding author: Aleksandar Lebl IRITEL a.d., 11080 Belgrade, Batajnički put 23, Serbia E-mail: lebl@iritel.com 380 A. LEBL, M. MILEUSNIĆ, D. MITIĆ, J. RADIVOJEVIĆ, V. MATIĆ The most often applied sensors for drone DILC are radar, RF signal detector, optical camera, thermal camera and acoustic detector. The benefits and drawbacks of each sensor type are emphasized in details in [2]. Drone DILC is usually performed using several sensor types among the mentioned ones. These selected sensors are then combined in one solution [3]-[6]. Among these sensor types radar, especially Frequency Modulated Continuous Wave (FMCW) radar, is probably the most often applied technique [7]-[17]. The main principles of FMCW radar realization are described in [18]-[20]. FMCW radar allows reliable classification of the detected drone based on the analysis of drone micro- Doppler signatures. Several typical drone construction and functional characteristics such as the number of its rotors, the number of blades in each rotor, rotor angular velocity and the length of blades may be only determined by FMCW radar on the base of drone micro- Doppler signature even in bad weather conditions. Among contributions in the domain of FMCW radar, micro-Doppler signatures for various drone types are presented and analyzed in [8]-[11], [15], [17]. Contribution [8] gives several drone micro-Doppler signature graphs in various flying mode phases (take- off, hovering, flying phase). In this aspect [8] is more complete than our paper, but it lacks explanation to make the connection between the graphs and the derived formulas for micro-Doppler signature calculation. The paper [9] presents a number of micro- Doppler signature graphs, but with addition of signals used for the communication between the drone and its operator, signals for drone video communication and so on. In [9] it is not possible to distinguish the spectrum behaviour as a consequence of drone flying from other frequency spectrum sources components. Contributions [10], [15] are interesting because they pave the way in the comparison of drone and birds micro-Doppler signatures, because drone and birds are often hard to distinguish due to their similar dimensions. The paper [21] contains very detailed theoretical and practical analysis of drone micro-Doppler signatures, but on the base of experiments performed for drones at the distance of only several meters from the radar. Drone micro-Doppler signature graphs are often analyzed applying artificial intelligence algorithms, as for example in [22]. Elements which have influence on the characteristics of drone micro-Doppler signature are briefly emphasized in the Section 2. The calculation method for drone micro-Doppler signature determination is described in detail in the Section 3. The calculation method is illustrated by a number of examples in the Section 4 when drone physical characteristics and position relative to FMCW radar are changed. The concluding comments are given in the Section 5. 2. DRONE PARTS CAUSING MICRO-DOPPLER EFFECT All drone moving parts may cause micro-Doppler effect detectable by FMCW radar. It is very important, because even a drone in hovering state will be detected by radar sensor. Drone micro-movable parts are its rotors. Each drone has a certain number of rotors, as presented in the Fig. 1. There are Nr=4 rotors in the example from Fig. 1. Drone micro-Doppler signature depends on this number of rotors. The second important factor which has influence on drone micro-Doppler signature is the number of blades (N) in each rotor. The blades 1 and 2 are designated in the Fig. 1. Verification of Calculation Method for Drone Micro-Doppler Signature Estimation 381 Two remaining blade characteristics which determine micro-Doppler signature are the blades length (L) and blades rotation speed (Ω). Drone micro-Doppler echo also depends on the drone (i.e. drone rotors) elevation angle (β) towards the radar level. This angle is determined by the drone height (h) and its distance from the radar (R0). h L Blade Blade Blade Blade L Fig. 1 Elements which have influence on the drone micro-Doppler signature 3. CALCULATION METHOD Method for drone micro-Doppler signature calculation may be explained referring to the Fig. 1. The main characteristic of FMCW radar is that it generates signal of variable frequency as a function of time. This frequency change is usually linear (sweep signal) and it is essentially important for FMCW radar detection principle of operation. The generated signal may be expressed by the equation [23] ( ) cos(2 ( ) ) c s t f b t t=  +   (1) where fc is the starting frequency of FMCW radar sweep signal and b is the slope of generated sweep signal. The generated signal is periodically repeated. The returned echo signal from rotor blades may be expressed by the equation from [7]: 1 0 0 0 1 0 4 ( ) ( ) exp ( sin ) sin ( ( )) exp( ( )) N lk k N k k k S t s t L j R z c t j t −  = − =   = =  −  +          −    (2) where it is 382 A. LEBL, M. MILEUSNIĆ, D. MITIĆ, J. RADIVOJEVIĆ, V. MATIĆ 0 4 2 ( ) cos cos ( =0,1,2,... -1). 2 k L k t t k N N     =     + +     (3) In these two equations: ▪ L is the length of each blade; ▪ N is the number of blades in each rotor; ▪ R0 is the distance between the radar and the drone rotor (approximately the same as between radar and drone); ▪ z0 is the drone height; ▪ β is the drone elevation angle in relation to radar; ▪ Ω is rotor angular rotation speed; ▪ φ0 is initial rotation angle; ▪ λ is FMCW signal wavelength. The magnitude of the rotor echo signal is 0 0 1 0 4 exp ( sin ) ( ) . sin ( ( )) exp( ( )) N k k k L j R z S t c t j t  − =    −  +       =    −  (4) The echo signal of all drone rotors is calculated according to the expression from [8]: 1 0 0 1 0 1 1 0 4 ( ) ( ) exp ( sin ) sin ( ( )) exp( ( )) r rN NN lk i i i i k i N ik ik k S t s t L j R z c t j t −  = = = − =   = =  −  +          −     (5) where Nr is the number of drone rotors and 0 4 2 ( ) cos cos ( =0,1,2,... -1). 2 ik i i i L k t t k N N     =       + +     (6) As for the case of only one rotor, the magnitude of the whole drone echo signal is, similar to the equation (4): 0 01 1 1 1 0 0 4 exp ( sin ) ( ) ( ) . sin ( ( )) exp( ( )) r r N i i iN N i lk N i k ik ik k L j R z S t s t c t j t − =  − = = =    −  +       = =    −     (7) The usual way to analyze drone micro-Doppler signatures is the application of drone spectrograms. Spectrograms present frequency spectrum of a signal as a function of time. They are obtained after calculation of Short-Time Fourier Transform (STFT) [24]: ( ( , )) exp( ) n n n m n n STFT S m S w j t  − =−  =   −  (8) or in the logarithmic division ( ) 20 log ( ( , )) n STFT dB STFT S m=   (9) The meaning of variables in (8) is: Verification of Calculation Method for Drone Micro-Doppler Signature Estimation 383 ▪ Sn – sequence of time samples of the signal whose spectrogram is calculated; ▪ wn – sequence of time samples of the selected window function; ▪ m – time index, i.e. time shift of the moment for which spectrogram is calculated; ▪ ω – frequency of the signal. Hanning window is most often selected for the calculation of STFT. The sequence of discretized Hanning windows function is expressed as [25]: 1 2 1 cos ( 0,1, 2... ) 2 n n w n N N   =  − =    (10) 4. DRONE SPECTROGRAMS Drone spectrograms obtained using equations (2) to (9) are presented in the figures 2 to 10. They are derived varying the mechanical and position characteristics of drones to analyze how the change of each parameter influences the spectrogram. The analysis is presented for the hovering drone which means that rotor blades are the only moving parts of the drone. The majority of spectrograms are presented for a single rotor and this corresponds to the class of drones in the shape of helicopter. This is the smaller in number class then the class in the shape of quadcopters (which have four rotors). The starting spectrogram is presented in the Fig. 2. It corresponds to the case that there is only one rotor with one blade. The blade rotation speed is Ωrot=30rotations/s and the blade length is L=0.24m. After these mechanical characteristics, the drone position in relation to radar is defined by its height h=30m and distance from radar R0=100m meaning that drone position elevation angle relative to radar is β=arc sin (0.3). The radar functional characteristics are operating frequency f=24GHz (operating wavelength 0.0125m) and sampling rate fstep=20kHz. Time interval for spectrogram presentation is 0.1s and waveform repetition rate during this time interval is 3 or 30 in 1s. It means that spectrogram appearance (time repetition rate) directly follows from the rotor rotation speed. 0 ,0 0 2 0 ,0 0 6 0 ,0 1 0 ,0 1 4 0 ,0 1 8 0 ,0 2 2 0 ,0 2 6 0 ,0 3 0 ,0 3 4 0 ,0 3 8 0 ,0 4 2 0 ,0 4 6 0 ,0 5 0 ,0 5 4 0 ,0 5 8 0 ,0 6 2 0 ,0 6 6 0 ,0 7 0 ,0 7 4 0 ,0 7 8 0 ,0 8 2 0 ,0 8 6 0 ,0 9 0 ,0 9 4 0 ,0 9 8 0 16 32 48 64 80 96 112 128 144 160 176 192 t (sec) f (Hz) 0-10 -10-0 -20--10 -30--20 -40--30 -50--40 -60--50 -70--60 -80--70 -90--80 -100--90 -110--100 -120--110 -130--120 -140--130 0 ,0 0 2 0 ,0 0 6 0 ,0 1 0 ,0 1 4 0 ,0 1 8 0 ,0 2 2 0 ,0 2 6 0 ,0 3 0 ,0 3 4 0 ,0 3 8 0 ,0 4 2 0 ,0 4 6 0 ,0 5 0 ,0 5 4 0 ,0 5 8 0 ,0 6 2 0 ,0 6 6 0 ,0 7 0 ,0 7 4 0 ,0 7 8 0 ,0 8 2 0 ,0 8 6 0 ,0 9 0 ,0 9 4 0 ,0 9 8 0 16 32 48 64 80 96 112 128 144 160 176 192 t (sec) f (Hz) 0-10 -10-0 -20--10 -30--20 -40--30 -50--40 -60--50 -70--60 -80--70 -90--80 -100--90 -110--100 -120--110 -130--120 -140--130 Fig. 2 Drone spectrogram for one rotor with one blade, the blade length L=0.24m, blade rotation speed Ωrot=30rotations/s, drone height h=30m, drone distance from radar R0=100m, FMCW radar operating frequency f=24GHz, digital sampling rate fstep=20kHz. 384 A. LEBL, M. MILEUSNIĆ, D. MITIĆ, J. RADIVOJEVIĆ, V. MATIĆ For our analysis in this paper it is important to notice the frequency at which signal echo falls below -40dB, i.e. where the spectrogram colour transfers from yellow to green. This frequency in the case of spectrogram from the Fig. 2 is 144Hz. Fig. 3 presents the drone spectrogram for the same parameters as in the Fig. 2 with the only difference that the blade rotation speed is Ωrot=20rotations/s. Two modifications are noticeable as a consequence of Ωrot change: the signal repetition rate has dropped from 3 to 2 during 0.1s and the frequency at which signal echo falls below -40dB is 96Hz. In both cases the parameter ratio is 2/3 as also the ratio of Ωrot values. This change of the important frequency bandwidth is important for our future analysis. Fig. 4 presents the drone spectrogram for the case that the blade length has been changed comparing to the Fig. 2. In this case the frequency at which signal echo falls below -40dB is a bit more than 72Hz. It means that in this case the ratio of important frequencies bandwidth has dropped in the ratio 1/2, as also the ratio of blades length. 0 ,0 0 2 0 ,0 0 6 0 ,0 1 0 ,0 1 4 0 ,0 1 8 0 ,0 2 2 0 ,0 2 6 0 ,0 3 0 ,0 3 4 0 ,0 3 8 0 ,0 4 2 0 ,0 4 6 0 ,0 5 0 ,0 5 4 0 ,0 5 8 0 ,0 6 2 0 ,0 6 6 0 ,0 7 0 ,0 7 4 0 ,0 7 8 0 ,0 8 2 0 ,0 8 6 0 ,0 9 0 ,0 9 4 0 ,0 9 8 0 16 32 48 64 80 96 112 128 144 160 176 192 t (sec) f (Hz) 0-10 -10-0 -20--10 -30--20 -40--30 -50--40 -60--50 -70--60 -80--70 -90--80 -100--90 -110--100 -120--110 -130--120 -140--130 0 ,0 0 2 0 ,0 0 6 0 ,0 1 0 ,0 1 4 0 ,0 1 8 0 ,0 2 2 0 ,0 2 6 0 ,0 3 0 ,0 3 4 0 ,0 3 8 0 ,0 4 2 0 ,0 4 6 0 ,0 5 0 ,0 5 4 0 ,0 5 8 0 ,0 6 2 0 ,0 6 6 0 ,0 7 0 ,0 7 4 0 ,0 7 8 0 ,0 8 2 0 ,0 8 6 0 ,0 9 0 ,0 9 4 0 ,0 9 8 0 16 32 48 64 80 96 112 128 144 160 176 192 t (sec) f (Hz) 0-10 -10-0 -20--10 -30--20 -40--30 -50--40 -60--50 -70--60 -80--70 -90--80 -100--90 -110--100 -120--110 -130--120 -140--130 Fig. 3 Drone spectrogram for one rotor with one blade, the blade length L=0.24m, blade rotation speed Ωrot=20rotations/s, drone height h=30m, drone distance from radar R0=100m, FMCW radar operating frequency f=24GHz, digital sampling rate fstep=20kHz. 0 ,0 0 2 0 ,0 0 6 0 ,0 1 0 ,0 1 4 0 ,0 1 8 0 ,0 2 2 0 ,0 2 6 0 ,0 3 0 ,0 3 4 0 ,0 3 8 0 ,0 4 2 0 ,0 4 6 0 ,0 5 0 ,0 5 4 0 ,0 5 8 0 ,0 6 2 0 ,0 6 6 0 ,0 7 0 ,0 7 4 0 ,0 7 8 0 ,0 8 2 0 ,0 8 6 0 ,0 9 0 ,0 9 4 0 ,0 9 8 0 16 32 48 64 80 96 112 128 144 160 176 192 t (sec) f (Hz) 0-10 -10-0 -20--10 -30--20 -40--30 -50--40 -60--50 -70--60 -80--70 -90--80 -100--90 -110--100 -120--110 -130--120 -140--130 0 ,0 0 2 0 ,0 0 6 0 ,0 1 0 ,0 1 4 0 ,0 1 8 0 ,0 2 2 0 ,0 2 6 0 ,0 3 0 ,0 3 4 0 ,0 3 8 0 ,0 4 2 0 ,0 4 6 0 ,0 5 0 ,0 5 4 0 ,0 5 8 0 ,0 6 2 0 ,0 6 6 0 ,0 7 0 ,0 7 4 0 ,0 7 8 0 ,0 8 2 0 ,0 8 6 0 ,0 9 0 ,0 9 4 0 ,0 9 8 0 16 32 48 64 80 96 112 128 144 160 176 192 t (sec) f (Hz) 0-10 -10-0 -20--10 -30--20 -40--30 -50--40 -60--50 -70--60 -80--70 -90--80 -100--90 -110--100 -120--110 -130--120 -140--130 Fig. 4 Drone spectrogram for one rotor with one blade, the blade length L=0.12m, blade rotation speed Ωrot=30rotations/s, drone height h=30m, drone distance from radar R0=100m, FMCW radar operating frequency f=24GHz, digital sampling rate fstep=20kHz. Verification of Calculation Method for Drone Micro-Doppler Signature Estimation 385 Fig. 5 presents the drone spectrogram for the case when its height has changed from h1=30m to h2=70m. It means that the ratio of elevation angle cosine functions has changed in the ratio 2 1 0 2 2 0 1 1.335 1 elev h R Q h R   −    = =   −    (11) The bandwidth of important frequencies has changed in approximately the same ratio: from 144Hz to about 109Hz for the limit of -40dB or, in other words, this ratio is 1.32. Fig. 6 presents the spectrogram for the more probable case that the rotor has two blades. The other parameters for this spectrogram are the same as in the Fig. 2. The important frequencies bandwidth remains 144Hz as in the Fig 2, but the repetition rate is twice as in the Fig. 2, or total 6 due to the increased number of blades. Highly similar spectrogram is obtained for the example of a rotor with one blade with two-fold rotation speed (Ωrot=60rotations/s) and half a blade length (L=0.12m) and special attention has to be paid to distinguish these two cases. The spectrogram for this second case is presented in the Fig. 7. 0 ,0 0 2 0 ,0 0 6 0 ,0 1 0 ,0 1 4 0 ,0 1 8 0 ,0 2 2 0 ,0 2 6 0 ,0 3 0 ,0 3 4 0 ,0 3 8 0 ,0 4 2 0 ,0 4 6 0 ,0 5 0 ,0 5 4 0 ,0 5 8 0 ,0 6 2 0 ,0 6 6 0 ,0 7 0 ,0 7 4 0 ,0 7 8 0 ,0 8 2 0 ,0 8 6 0 ,0 9 0 ,0 9 4 0 ,0 9 8 0 16 32 48 64 80 96 112 128 144 160 176 192 t (sec) f (Hz) 0-10 -10-0 -20--10 -30--20 -40--30 -50--40 -60--50 -70--60 -80--70 -90--80 -100--90 -110--100 -120--110 -130--120 -140--130 0 ,0 0 2 0 ,0 0 6 0 ,0 1 0 ,0 1 4 0 ,0 1 8 0 ,0 2 2 0 ,0 2 6 0 ,0 3 0 ,0 3 4 0 ,0 3 8 0 ,0 4 2 0 ,0 4 6 0 ,0 5 0 ,0 5 4 0 ,0 5 8 0 ,0 6 2 0 ,0 6 6 0 ,0 7 0 ,0 7 4 0 ,0 7 8 0 ,0 8 2 0 ,0 8 6 0 ,0 9 0 ,0 9 4 0 ,0 9 8 0 16 32 48 64 80 96 112 128 144 160 176 192 t (sec) f (Hz) 0-10 -10-0 -20--10 -30--20 -40--30 -50--40 -60--50 -70--60 -80--70 -90--80 -100--90 -110--100 -120--110 -130--120 -140--130 Fig. 5 Drone spectrogram for one rotor with one blade, the blade length L=0.24m, blade rotation speed Ωrot=30rotations/s, drone height h=70m, drone distance from radar R0=100m, FMCW radar operating frequency f=24GHz, digital sampling rate fstep=20kHz. 386 A. LEBL, M. MILEUSNIĆ, D. MITIĆ, J. RADIVOJEVIĆ, V. MATIĆ 0 ,0 0 2 0 ,0 0 6 0 ,0 1 0 ,0 1 4 0 ,0 1 8 0 ,0 2 2 0 ,0 2 6 0 ,0 3 0 ,0 3 4 0 ,0 3 8 0 ,0 4 2 0 ,0 4 6 0 ,0 5 0 ,0 5 4 0 ,0 5 8 0 ,0 6 2 0 ,0 6 6 0 ,0 7 0 ,0 7 4 0 ,0 7 8 0 ,0 8 2 0 ,0 8 6 0 ,0 9 0 ,0 9 4 0 ,0 9 8 0 16 32 48 64 80 96 112 128 144 160 176 192 t (sec) f (Hz) 0-10 -10-0 -20--10 -30--20 -40--30 -50--40 -60--50 -70--60 -80--70 -90--80 -100--90 -110--100 -120--110 -130--120 -140--130 0 ,0 0 2 0 ,0 0 6 0 ,0 1 0 ,0 1 4 0 ,0 1 8 0 ,0 2 2 0 ,0 2 6 0 ,0 3 0 ,0 3 4 0 ,0 3 8 0 ,0 4 2 0 ,0 4 6 0 ,0 5 0 ,0 5 4 0 ,0 5 8 0 ,0 6 2 0 ,0 6 6 0 ,0 7 0 ,0 7 4 0 ,0 7 8 0 ,0 8 2 0 ,0 8 6 0 ,0 9 0 ,0 9 4 0 ,0 9 8 0 16 32 48 64 80 96 112 128 144 160 176 192 t (sec) f (Hz) 0-10 -10-0 -20--10 -30--20 -40--30 -50--40 -60--50 -70--60 -80--70 -90--80 -100--90 -110--100 -120--110 -130--120 -140--130 Fig. 6 Drone spectrogram for one rotor with two blades, the blade length L=0.24m, blade rotation speed Ωrot=30rotations/s, drone height h=30m, drone distance from radar R0=100m, FMCW radar operating frequency f=24GHz, digital sampling rate fstep=20kHz. 0 ,0 0 2 0 ,0 0 6 0 ,0 1 0 ,0 1 4 0 ,0 1 8 0 ,0 2 2 0 ,0 2 6 0 ,0 3 0 ,0 3 4 0 ,0 3 8 0 ,0 4 2 0 ,0 4 6 0 ,0 5 0 ,0 5 4 0 ,0 5 8 0 ,0 6 2 0 ,0 6 6 0 ,0 7 0 ,0 7 4 0 ,0 7 8 0 ,0 8 2 0 ,0 8 6 0 ,0 9 0 ,0 9 4 0 ,0 9 8 0 16 32 48 64 80 96 112 128 144 160 176 192 t (sec) f (Hz) 0-10 -10-0 -20--10 -30--20 -40--30 -50--40 -60--50 -70--60 -80--70 -90--80 -100--90 -110--100 -120--110 -130--120 -140--130 0 ,0 0 2 0 ,0 0 6 0 ,0 1 0 ,0 1 4 0 ,0 1 8 0 ,0 2 2 0 ,0 2 6 0 ,0 3 0 ,0 3 4 0 ,0 3 8 0 ,0 4 2 0 ,0 4 6 0 ,0 5 0 ,0 5 4 0 ,0 5 8 0 ,0 6 2 0 ,0 6 6 0 ,0 7 0 ,0 7 4 0 ,0 7 8 0 ,0 8 2 0 ,0 8 6 0 ,0 9 0 ,0 9 4 0 ,0 9 8 0 16 32 48 64 80 96 112 128 144 160 176 192 t (sec) f (Hz) 0-10 -10-0 -20--10 -30--20 -40--30 -50--40 -60--50 -70--60 -80--70 -90--80 -100--90 -110--100 -120--110 -130--120 -140--130 Fig. 7 Drone spectrogram for one rotor with one blade, the blade length L=0.12m, blade rotation speed Ωrot=60rotations/s, drone height h=30m, drone distance from radar R0=100m, FMCW radar operating frequency f=24GHz, digital sampling rate fstep=20kHz. The echo signal at the frequency 0Hz may be used to distinguish whether it is considered the case according to the Fig. 6 or the Fig. 7. Echo signal amplitude oscillations are significantly greater when rotation speed is lower, as is illustrated by the characteristics presented in the Fig. 8 and the Fig. 9. This peak-to-peak amplitude of the oscillations is even about 17dB when there are two blades of 0.24m length and their rotation speed is 30 rotations/s (Fig. 8) comparing to only about 2.5dB when there is one blade of 0.12m length moving at Ωrot=60rotations/s (Fig. 9). This presentation of echo Verification of Calculation Method for Drone Micro-Doppler Signature Estimation 387 signal at the frequency 0Hz for spectrograms more reliable distinguishing in some cases is, as for our knowledge, the paper original contribution. -30 -25 -20 -15 -10 -5 0 0,002 0,012 0,022 0,032 0,042 0,052 0,062 0,072 0,082 0,092 t [s] a [ d B ] Fig. 8 Echo at the frequency 0Hz for the case of one rotor with two blades, the blade length L=0.24m, blade rotation speed Ωrot=30rotations/s, drone height h=30m, drone distance from radar R0=100m, FMCW radar operating frequency f=24GHz, digital sampling rate fstep=20kHz. -30 -25 -20 -15 -10 -5 0 0,002 0,012 0,022 0,032 0,042 0,052 0,062 0,072 0,082 0,092 t [s] a [ d B ] Fig. 9 Echo at the frequency 0Hz for the case of one rotor with one blade, the blade length L=0.12m, blade rotation speed Ωrot=60rotations/s, drone height h=30m, drone distance from radar R0=100m, FMCW radar operating frequency f=24GHz, digital sampling rate fstep=20kHz. The typical drone construction is with 4 rotors and each rotor with two blades. The spectrogram for such a construction is presented in the Fig. 10. The consequence of more rotors and blades existence is that echo signal periodicity is less obvious and that limit value of important echo frequencies is practically constant as a function of time. 388 A. LEBL, M. MILEUSNIĆ, D. MITIĆ, J. RADIVOJEVIĆ, V. MATIĆ 0 ,0 0 2 0 ,0 0 6 0 ,0 1 0 ,0 1 4 0 ,0 1 8 0 ,0 2 2 0 ,0 2 6 0 ,0 3 0 ,0 3 4 0 ,0 3 8 0 ,0 4 2 0 ,0 4 6 0 ,0 5 0 ,0 5 4 0 ,0 5 8 0 ,0 6 2 0 ,0 6 6 0 ,0 7 0 ,0 7 4 0 ,0 7 8 0 ,0 8 2 0 ,0 8 6 0 ,0 9 0 ,0 9 4 0 ,0 9 8 0 16 32 48 64 80 96 112 128 144 160 176 192 t (sec) f (Hz) 0-10 -10-0 -20--10 -30--20 -40--30 -50--40 -60--50 -70--60 -80--70 -90--80 -100--90 -110--100 -120--110 -130--120 -140--130 0 ,0 0 2 0 ,0 0 6 0 ,0 1 0 ,0 1 4 0 ,0 1 8 0 ,0 2 2 0 ,0 2 6 0 ,0 3 0 ,0 3 4 0 ,0 3 8 0 ,0 4 2 0 ,0 4 6 0 ,0 5 0 ,0 5 4 0 ,0 5 8 0 ,0 6 2 0 ,0 6 6 0 ,0 7 0 ,0 7 4 0 ,0 7 8 0 ,0 8 2 0 ,0 8 6 0 ,0 9 0 ,0 9 4 0 ,0 9 8 0 16 32 48 64 80 96 112 128 144 160 176 192 t (sec) f (Hz) 0-10 -10-0 -20--10 -30--20 -40--30 -50--40 -60--50 -70--60 -80--70 -90--80 -100--90 -110--100 -120--110 -130--120 -140--130 Fig. 10 Drone spectrogram for four rotors with two blades, the blade length L=0.24m, blade rotation speed Ωrot=30rotations/s, drone height h=30m, drone distance from radar R0=100m, FMCW radar operating frequency f=24GHz, digital sampling rate fstep=20kHz. The graphs in the figures 2-7 may be compared to the selected graph from [8] which corresponds to the micro-Doppler signature of rotors obtained by practical recording. The great similarity is obvious with the exception that the graph in [8] is presented for positive and negative frequencies and the echo signal is symmetrical about the frequency 0Hz. This graph from [8] is presented in the Fig. 11. It has periodicity – the number of periodical changes is 18 during 1s. According to this characteristic, the graph is most similar to the graph in the Fig. 3. The frequency where the signal echo rapidly decreases is about 100Hz. Let us further suppose that we could conclude by some other technique what is the drone elevation angle, i.e. what is cosine of elevation angle. The final element to determine is now the length of rotor blade/blades (L). Under the assumption that elevation angle is the same as in the Fig. 3, we obtain L=0.25m. But, if the drone is situated approximately vertically above the FMCW radar (i.e. elevation angle tends to 90o) and the spectrogram is without significant changes, the corresponding L quickly grows. The graph in the Fig. 10 is similar to the graph from [8] which corresponds to the drone in the hovering state. This graph from [8] is presented in the Fig. 12. There is no obvious periodicity in the recorded characteristic. Such a graph is the clear sign that there is a higher number of rotors probably with more than one blade. The summary of conditions for spectrogram characteristic calculation in the figures 2-10 is presented in the Table 1. The main specificities to describe the obtained spectrograms for each combination of conditions (i.e. each figure) are also presented in the Table 1. Verification of Calculation Method for Drone Micro-Doppler Signature Estimation 389 Fig. 11 Practical rotor micro-Doppler record [8] Fig. 12 Practical drone micro-Doppler record [8] Table 1 Summary of figure characteristics Figure Conditions for spectrogram calculation Output spectrogram description 2 1 rotor, 1 blade, L=0.24m, Ωrot=30/s, h=30m, R0=100m, f=24GHz, fstep=20kHz Waveform repetition rate 30/s Attenuation 40dB at 144Hz 3 Figure 2 with Ωrot=20/s Waveform repetition rate 20/s Attenuation 40dB at 96Hz 4 Figure 2 with L=0.12/s Waveform repetition rate 30/s Attenuation 40dB at 72Hz 5 Figure 2 with h=70m (cosine of elevation angle higher 1.335 times) Waveform repetition rate 30/s Attenuation 40dB at 109Hz 6 Figure 2 with two blades Waveform repetition rate 60/s Attenuation 40dB at 144Hz 7 Figure 2 with L=0.12/s and Ωrot=60/s Waveform repetition rate 60/s Attenuation 40dB at 144Hz 8 Figure 2 with two blades Amplitude oscillations peak-to-peak 17dB at 0Hz 9 Figure 2 with L=0.12/s and Ωrot=60/s Amplitude oscillations peak-to-peak 2.5dB at 0Hz 10 Figure 2 with four rotors and two blades Echo frequencies constant in time, signal periodicity less obvious 390 A. LEBL, M. MILEUSNIĆ, D. MITIĆ, J. RADIVOJEVIĆ, V. MATIĆ 5. CONCLUSIONS Calculation method for drone micro-Doppler signature determination is presented in this paper. The influence of various drone parameters (number of rotors, number of blades forming a rotor, blades rotation rate, blades length) on spectrogram shape is analyzed. Special attention is devoted to the way how it is possible to distinguish some combinations of drone characteristics which give very similar spectrograms. All results are presented for the FMCW radar which operates on the frequency of 24GHz. The method and the results from the paper may be used in the case that measurement results are not available. The results of calculation are compared to the similar examples from measurements and similarity of the results from these two groups is verified by several practical examples. The results from this paper are related only to the hovering drone. Our plan for the future investigation is to try to develop calculation method for the drones in other flying modes (flying, take-off, etc). Multi-Doppler spectrograms are applicable for drone detection, identification and classification by artificial intelligence algorithms. Our other development direction plan is to implement calculated spectrograms for training neural networks in the first phase of such networks construction when numerous practical records of various drone types are still not available. REFERENCES [1] V. Matić, V. Kosjer, A. Lebl, B. Pavić and J. Radivojević, "Methods for Drone Detection and Jamming", In Proceedings of the 10th International Conference on Information Society and Technology (ICIST). Kopaonik, 2020, pp.16–21. [2] N. Eriksson, Conceptual study of a future drone detection system Countering a threat posed by a disruptive technology. Master thesis in Product Development, Chalmers University of Technology, Goethenburg, Sweden, 2018. [3] Advanced protection systems, Ctrl+sky drone detection and neutralization system, 2017, http://apsystems.tech/wp-content/uploads/2018/01/aps_broszura_web.pdf. [4] Droneshield, "Product Information", 2018. [5] H. Liu, F. Qu, Y. Liu, W. Zhao and Y. Chen, "A drone detection with aircraft classification based on a camera array", In Proceedings of the 2018 IOP Conference Series: Materials Science and Engineering, vol. 322, p. 052005. 2018, pp. 1–7. [6] X. Shi, C. Yang, C. Liang, Z. Shi and J. Chen, "Anti-Drone System with Multiple Surveillance Technologies: Architecture, Implementation, and Challenges", IEEE Commun. Magaz., vol. 56, no. 4, pp. 68–74, 2018. [7] V. C. Chen, The Micro-Doppler Effect in Radar. Artech House, Second Edition, 2019, ISBN: 978-1-63081- 546-2. [8] C. Zhao, G. Luo, Y. Wang, C. Chen and Z. Wu, "UAV Recognition Based on Micro-Doppler Dynamic Attribute-Guided Augmentation Algorithm", Remote Sensing, vol. 13, no. 6, p. 1205, pp. 1–17, 2021. [9] T. Šević, V. Joksimović, I. Pokrajac, R. Brusin, B. Sazdić-Jotić and D. Obradović, "Interception and Detection of Drones Using RF-based Dataset of Drones", Sci. Tech. Rev., vol. 70, no. 2, pp. 29–34, 2020. [10] S. Rahman and D. A. Robertson, "Radar micro-Doppler signatures of drones and birds at K-band and W- band", Sci. Rep., vol. 8, pp. 1–11, 2018. [11] Y. Cai, O. Krasnov and A. Yarovoy, "Simulation of Radar Micro-Doppler Patterns for Multi-propeller Drones", In Proceedings of the International Radar Conference (Radar-2019), Toulon, 2019, pp.1–5. [12] W. Wang, J. Du and J. Gao, "Multi-Target Detection Method Based on Variable Carrier Frequency Chirp Sequence", Sensors, vol. 18, p. 3386, pp. 1–12, 2018. [13] A. Coluccia, G. Parisi and A. Fascista, "Detection and Classification of Multirotor Drones in Radar Sensor Networks: A Review", Sensors, vol. 20, p. 4172, pp. 1–22, 2020. http://apsystems.tech/wp-content/uploads/2018/01/aps_broszura_web.pdf Verification of Calculation Method for Drone Micro-Doppler Signature Estimation 391 [14] M. Daković, M. Brajović, T. Thayaparan and Lj. Stanković, "An algorithm for micro-Doppler period estimation", In Proceedings of the 20th Telecommunications Forum (TELFOR), Belgrade, 2012, pp. 851–854. [15] P. Molchanov, Radar Target Classification by Micro-Doppler Contributions. Thesis for the degree of Doctor of Science in Technology, Publication 1255, Tampere University of Technology, Finland, October 2014, ISSN 1459-2045. [16] E. Hyun, Y.-S. Jin and J.-H. Lee, "Design and Implementation of 24 GHz Multichannel FMCW Surveillance Radar with a Software-Reconfigurable Baseband", J. Sensors, vol. 2017, p. 3148237, pp. 1–11, 2017. [17] B. Karlsson, Modeling multicopter radar return. Master’s thesis in Applied Physics, Chalmers University of Technology, Department of Electrical Engineering, Gothenburg, Sweden, 2017. [18] V. M. Milovanović, “On Fundamental Operating Principles and Range-Doppler Estimation in Monolithic Frequency-Modulated Continuous-Wave Radar Sensors", FU Elec. Energ., vol. 31, no. 4, pp. 547–570, 2018. [19] C. Iovescu and S. Rao, The fundamental of millimeter wave radar sensors. Texas Instruments, 2020. [20] J. Zhu, Low-cost, software defined FMCW radar for observations of drones. Master thesis, University of Oklahoma, Graduate College, 2017. [21] M. Passafiume, N. Rojhani, G. Collodi and A. Cidronali, "Modeling Small UAV Micro-Doppler Signature Using Milimeter-Wave FMCW Radar", Electronics , vol. 10, no. 6, pp. 1–16, 2021. [22] J. Park, J.-S. Park and S.-O. Park, "Small Drone Classification with Light CNN and New Micro-Doppler Signature Extraction Method Based on A-SPC Technique", https://arxiv.org/abs/2009.14422, pp.1–5, 2020. [23] T. Tang and C. Wu, Design of new Frequency Modulated Continuous Wave (FMCW) target tracking radar with digital beamforming tracking. Defense Research and Development Canada, Scientific Report DRDC- RDDC-2019-R175, 2019. [24] M. Ahmadizadeh, An Introduction to Short-Time Fourier Transform (STFT). Sharif University of Technology, Department of Civil Engineering, 2014. [25] H. A. Gaberson, "A Comprehensive Windows Tutorial", Sound and Vibration, Instrumentation Reference Issue, pp. 14–23, 2006. https://arxiv.org/abs/2009.14422