[ADD] Basic FM IoO simulator, range-Doppler matrix slice plot

README.md

@@ -1,10 +1,11 @@
 # py A          P         R      i  L
+
 # Advanced   Passive   Radar     Library
 
 pyAPRiL is a python based signal processing library which implements passive radar signal processing algorithms. All the algorithms are tested and verified through real field measurement data and simulations. The corresponding references are highlited where applicable. Algorithms are researched and developed on real life DVB-T and FM
 based passive radar systems.
 
-##### The package is organized as follows:
+### The package is organized as follows:
 
 * pyAPRiL: Main package directory
 	* **channelPreparation**: This file contains wrapper functions for the different algorithms that aims to prepare the  reference and the sureveillance signals for the detection stage. Such as reference signal regeneration and clutter cancellation techniques.
@@ -39,31 +40,32 @@ based passive radar systems.
 	    * P: Peak clutter reduction
 	    * D: Dynamic range compression
 	* **targetParameterCalculator**: Calculates the expected parameters of an observed target (bistatic range, Doppler, ...) from geodetic data.
-	* **RDTools**: range-Doppler matrix image export tool
+	* **RDTools**: range-Doppler matrix image export and plotting tools
+	* **sim**: IoO simulator (FM)
 	* **docs**: Contains Ipython notebook files with demonstrations. (Currently not available)
 	* **testing**: Contains demonstration functions.
 
-##### Installing from Python Package Index:
+### Installing from Python Package Index:
 ```python
 pip install pyapril
 ```
 
-##### Version history
+### Version history
 * 1.0.0 Inital version
 * 1.1.0 Automatic detection added (CA-CFAR)
 * 1.2.0 Add temporal clutter filters (ECA, ECA-B, ECA-S)
 * 1.3.0 Clutter filter and detection performance evaluation
 * 1.4.0 Target DoA estimation
 * 1.5.0 Surveillance channel beamforming
 * 1.6.0 Target parameter calculator
-
-##### Acknowledgements
+* 1.7.0 FM IoO simulator
+### Acknowledgements
 This work was supported by the Microwave Remote Sensing Laboratory of BME ([Radarlab](http://radarlab.mht.bme.hu)). A special thanks of mine goes to the  RTL-SDR site ([RTL-SDR](https://www.rtl-sdr.com/)) who helped this project becoming mature with the development and testing of the Kerberos SDR.
 
 
 For further information on passive radars check: [tamaspeto.com](https://tamaspeto.com)
 *Tamás Pető*
-*2017-2020, Hungary*
+*2017-2021, Hungary*
 
 
 

pyapril/RDTools.py

@@ -273,3 +273,185 @@ def plot_rd_matrix(rd_matrix,
 
 
     return fig
+
+def plot_Doppler_slice(rd_matrix, bistat_range, **kwargs):
+    """
+        Description:
+        ------------
+        Displays the requested range slice of a given range-Doppler map.
+        In case the sampling frequency is specified through the kwargs parameter ('fs')
+        the requested range parameter is interpreted as bistatic range in [m], otherwise
+        it is interpreted as [bin]
+        
+        When a Plotly combatible figure object is passed through the 'fig' keyword,
+        the function will plot the extracted slice onto this figure, thus enabling
+        multiple slices to plot on the same figure object.
+        
+        
+        Parameters:
+        -----------
+        :param: rd_matrix              : range-Doppler matrix to be exported        
+        :param: bistat_range           : range index, either [bin] or [m]
+            
+        :type: rd_matrix              : R x D complex numpy array
+        :type: bistat_range           : float
+        
+         **kwargs
+        Additional display option can be specified throught the **kwargs interface
+        Valid keys are the followings:
+        
+        
+        :key:           fs: sampling frequency of the processed signal 
+        :key:  max_Doppler: Maximum Doppler frequency
+        :key:          fig: Figure objet to plot on
+
+        :type:          fs: float
+        :type: max_Doppler: float
+        :type:         fig: Plotly figure object
+       
+        
+        Return values:
+        --------------
+        :return: fig: Generated Doppler slice figure
+        :rtype:  fig: Plotly compatibile Figure object
+    
+    """    
+
+    """
+    --------------------------------
+               Parameters
+    --------------------------------
+    """
+
+    fs          = kwargs.get('fs')     
+    max_Doppler = kwargs.get('max_Doppler')    
+    fig         = kwargs.get('fig')
+    rd_matrix = 10 * np.log10(np.abs(rd_matrix) ** 2)    
+
+
+    """
+    --------------------------------
+            Generate Figure
+    --------------------------------
+    """    
+    if bistat_range < 0:
+        print("ERROR: Bistatic range should be a positive number{:.1f}".format(bistat_range))
+        return None
+
+    if fs is not None:
+        d =  3*10**8/fs
+        if bistat_range > (np.size(rd_matrix,1)-1)*d:
+            print("ERROR: Bistatic range is out of range. {:.1f} > Max bistatic range: {:.1f} m".format(bistat_range, np.size(rd_matrix,1)*d))
+            return None
+
+        bistat_range=np.argmin(abs((np.arange(np.size(rd_matrix,1))*d-bistat_range)))
+
+    if bistat_range > np.size(rd_matrix,1)-1:
+        print("ERROR: Bistatic range is out of range. {:.1f}, Range size: {:.1f} bin".format(bistat_range, np.size(rd_matrix,1)))
+        return None
+
+    # Prepare scales
+    doppler_scale = np.arange(np.size(rd_matrix,0), dtype=float)
+    name          = "Bistatic range bin: {:d}".format(bistat_range)
+    x_axis_title = "Doppler frequency [bin]"
+    if max_Doppler is not None:
+        doppler_scale = np.linspace(-max_Doppler,max_Doppler, np.size(rd_matrix,0))
+        x_axis_title = "Doppler frequency [Hz]"
+
+
+    # Prepare figure object
+    if fig is None:        
+        fig = go.Figure()      
+
+    fig.add_trace(go.Scatter(x=doppler_scale, y=rd_matrix[:,bistat_range], name=name))
+    fig.update_xaxes(title_text=x_axis_title)
+    fig.update_yaxes(title_text="Amplitude [dB]")
+
+    return fig
+
+def plot_range_slice(rd_matrix, Doppler_freq, **kwargs):
+    """
+        Description:
+        ------------
+        Displays the requested range slice of a given range-Doppler map.
+        In case the Doppler resolutionis specified through the kwargs parameter ('fD_res')
+        the requested Doppler slice parameter is interpreted as Hz, otherwise
+        it is interpreted as [bin]
+        
+        When a Plotly combatible figure object is passed through the 'fig' keyword,
+        the function will plot the extracted slice onto this figure, thus enabling
+        multiple slices to plot on the same figure object.
+        
+        Parameters:
+        -----------
+        :param: rd_matrix              : range-Doppler matrix to be exported        
+        :param: Doppler_freq           : Doppler frequency index, either [bin] or [Hz]
+            
+        :type: rd_matrix              : R x D complex numpy array
+        :type: Doppler_freq           : float
+        
+         **kwargs
+        Additional display option can be specified throught the **kwargs interface
+        Valid keys are the followings:
+        
+        
+        :key:   fs: sampling frequency of the processed signal   
+        :key:  fig: Figure objet to plot on
+        :type:  fs: float
+        :type: fig: Plotly figure object
+        
+        
+        Return values:
+        --------------
+        :return: fig: Generated range slice figure
+        :rtype:  fig: Plotly compatibile Figure object
+    
+    """    
+
+    """
+    --------------------------------
+               Parameters
+    --------------------------------
+    """
+
+    fs          = kwargs.get('fs')
+    max_Doppler = kwargs.get('max_Doppler')    
+    fig         = kwargs.get('fig')    
+    rd_matrix = 10 * np.log10(np.abs(rd_matrix) ** 2)    
+
+    """
+    --------------------------------
+            Generate Figure
+    --------------------------------
+    """
+    if max_Doppler is not None:
+        if abs(Doppler_freq) > max_Doppler:
+            print("ERROR: Doppler frequency is out of range. {:.1f} > Max Doppler: {:.1f} Hz".format(Doppler_freq, max_Doppler))
+            return None
+
+        Doppler_scale = np.linspace(-max_Doppler,max_Doppler, np.size(rd_matrix,0))
+        Doppler_freq=np.argmin(abs((Doppler_scale-Doppler_freq)))
+
+    if Doppler_freq > np.size(rd_matrix,0)-1:
+        print("ERROR: Doppler frequency is out of range. {:d}, Doppler size: {:d} bin".format(Doppler_freq, np.size(rd_matrix,0)))
+        return None
+
+    if Doppler_freq < 0:
+        print("ERROR: Doppler frequency is out of range. {:d}, Doppler size: {:d} bin".format(Doppler_freq, np.size(rd_matrix,0)))
+        return None
+
+    bistat_range_scale = np.arange(np.size(rd_matrix,1), dtype=float)
+    name = "Doppler bin: {:d}".format(Doppler_freq)
+    x_axis_title = "Bistatic range [bin]"
+    if fs is not None:
+        bistat_range_scale *= (3*10**8/fs)/10**3
+        x_axis_title = "Bistatic range [km]"
+
+    # Prepare figure object
+    if fig is None:        
+        fig = go.Figure()      
+    fig.add_trace(go.Scatter(x=bistat_range_scale, y=rd_matrix[Doppler_freq,:], name=name))
+    fig.update_xaxes(title_text=x_axis_title)
+    fig.update_yaxes(title_text="Amplitude [dB]")
+
+    return fig

pyapril/sim/FMSimulator.py

@@ -0,0 +1,97 @@
+import numpy as np
+import math
+import sys
+from pyapril.tools import resample
+
+from pydub import AudioSegment
+from pydub.utils import mediainfo
+
+"""
+                          Python based Advanced Passive Radar Library (pyAPRiL)
+
+                                            FM IoO Simulator
+"""
+
+def generate_fm_signal(s_mod, fs): 
+    """
+        Generates FM signal using the give modulation signal
+        
+        Parameters:
+        -----------
+        :param: s_mod: Modulating signal
+        :param:    fs: Sampling frequency of the modulating signal [Hz]
+        
+        :type: s_mod: One dimensional numpy array
+        :type:    fs: float
+        
+        Retrun values:
+        --------------
+        :return: s : Generated FM signal
+        :return: fs: Sampling frequency of the modulated signal
+        
+        :rtype:   s: One dimensional numpy array
+        :rtype:  fs: float 
+    
+    """
+
+    # Internal processing parameters
+    fd             = 75 # frequency deviation [kHz] -> Same as in FM broadcast
+
+    # Generate FM signal
+    s_mod = s_mod/np.max(np.abs(s_mod)) # Normalize
+    k_fm = fd*10**3/np.max(s_mod)    
+    s = np.sin(2*np.pi*k_fm/fs*np.cumsum(s_mod)) # Modulate
+    return s, fs
+
+def generate_fm_signal_from_sound(sound_fname, T_sim, offset): 
+    """
+    Generated FM signal from sound file
+    
+    The offset parameter is a float number between 0 and 1. It
+    specifyies the start position of the time window(T_sim) 
+    that will be used from the sound file for the FM signal 
+    preparation.
+    
+    Parameters:
+    ----------
+    :param: sound_fname: Name of sound sound file to be imported 
+    :param:       T_sim: Duration of the reqested simulated signal [s]
+    :param:      offset: Offset position of the processed window
+    
+    :type: sound_fname: string
+    :type:      T_sim : float
+    :type:      offset: float (0..1)
+    
+    Retrun values:
+    --------------
+    :return: s : Generated FM signal
+    :return: fs: Sampling frequency of the modulated signal
+    
+    :rtype:   s: One dimensional numpy array
+    :rtype:  fs: float 
+    
+    """
+    # Internal processing parameters
+    resample_ratio = 5 
+    fir_tap_size   = 10 # Resample filter tap size
+    fd             = 75 # frequency deviation [kHz] -> Same as in FM broadcast
+
+    # Import sound file
+    sound = AudioSegment.from_mp3(sound_fname)
+    sound_samples = np.array(sound.get_array_of_samples())
+
+    info     = mediainfo(sound_fname)
+    fs = int(info['sample_rate'])
+
+    # Resample    
+    offset = int(offset*np.size(sound_samples))
+    N_raw = math.ceil(T_sim/(1/fs)) # Number of samples needed from the modulating signal
+    sound_samples = sound_samples[offset:offset+N_raw] # Cut out useful portion
+    s_mod = resample.resample_n_filt(resample_ratio,1,fir_tap_size,sound_samples)
+    fs = resample_ratio * fs 
+
+    # Generate FM signal
+    s_mod = s_mod/np.max(np.abs(s_mod)) # Normalize
+    k_fm = fd*10**3/np.max(s_mod)
+    s = np.sin(2*np.pi*k_fm/fs*np.cumsum(s_mod)) # Modulate
+    return s, fs

pyapril/sim/__init__.py

@@ -0,0 +1 @@
+name = "pyapril"

pyapril/tools/__init__.py

@@ -0,0 +1 @@
+name = "pyapril"