Merge pull request #12 from AUAS-Pulsar/filterbankplot

Filterbankplot

Author: TVH7

README.md

@@ -7,11 +7,19 @@
 
 Creating a free and open source framework that contains the generic algorithms and file handling for astronomical data sets. This framework will be modular. Similar to OpenCV, wherein specific modules can be added and disabled depended on the needs of a project. This framework will be implemented in Python and C++.
 
-## Requirements
+
+
+
+## Installation
+
+### Requirements
 
     * numpy
     * python 3.6
 
+## Documentation
+You can find the documentation [here.](docs/README.md)
+
 # Table of contents
 1. [Getting started](docs/gettingstarted.md)
     1. [Example filterbank files](docs/gettingstarted.md#11-example-filterbank-files)

data_import.py

@@ -0,0 +1,25 @@
+import filterbank.filterbank as fb
+
+# Download the files in the README.md
+
+def _8bit(freq_start, freq_stop):
+    fb1 = fb.Filterbank(filename = "../data/pspm8.fil")
+
+    test1 = fb1.select_data(freq_start, freq_stop)
+
+    return test1
+
+
+def _16bit(freq_start, freq_stop):
+    fb2 = fb.Filterbank(filename = '../data/pspm16.fil')
+
+    test2 = fb2.select_data(freq_start, freq_stop)
+
+    return test2
+
+def _32bit(freq_start, freq_stop):
+    fb3 = fb.Filterbank(filename = '../data/pspm32.fil')
+
+    test3 = fb3.select_data(freq_start, freq_stop)
+
+    return test3

docs/README.md

@@ -0,0 +1,7 @@
+# Documentation
+## Filterbank
+[Filterbank docs](filterbank.md)
+## Fourier transforms
+[Fourier transforms docs](fourier.md)
+## Plotting
+[Plotting docs](plots.md)

docs/fourier.md

@@ -0,0 +1,40 @@
+# Fourier Transforms
+## DFT
+The `fourier.dft_slow` function is a plain implementation of discrete Fourier transformation. 
+
+The `fourier.fft_vectorized` function depends on this function.
+
+### Parameters
+
+| Parameter | Description |
+| --- | --- |
+| input_data | Array containing the values to be transformed. |
+
+Returns an array containing the transformed values.
+
+### Example usage: 
+```python
+>>> from Asteria import fourier
+>>> fourier.dft_slow([1,2,3,4])
+array([10.+0.00000000e+00j, -2.+2.00000000e+00j, -2.-9.79717439e-16j, -2.-2.00000000e+00j])
+```
+
+## FFT
+The `fourier.fft_vectorized` function is a vectorized, non-recursive version of the Cooley-Tukey FFT
+
+Gives the same result as `fourier.dft_slow` but is many times faster. 
+
+### Parameters:
+
+| Parameter | Description |
+| --- | --- |
+| input_data | Array containing the values to be transformed. |
+
+Returns an array containing the transformed values.
+
+### Example usage: 
+```python
+>>> from Asteria import fourier
+>>> fourier.fft_vectorized([1,2,3,4])
+array([10.+0.00000000e+00j, -2.+2.00000000e+00j, -2.-9.79717439e-16j, -2.-2.00000000e+00j])
+```

docs/plots.md

@@ -0,0 +1,60 @@
+# Plotting diagrams
+## PSD
+The `plot.psd` function can be used to generate the data for a Power Spectral Density plot. 
+
+### Parameters
+
+| Parameter | Description |
+| --- | --- |
+| samples | 1-D array or sequence. Array or sequence containing the data to be plotted. |
+| nfft | Integer, optional. The number of bins to be used. Defaults to 256. |
+| sample_rate | Integer, optional. The sample rate of the input data in `samples`. Defaults to 2.|
+| window | Callable, optional. The window function to be used. Defaults to `plot.window_hanning`. |
+| sides | {'default', 'onesided', 'twosided'}. Specifies which sides of the spectrum to return. Default gives the default behavior, which returns one-sided for real data and both for complex data. 'onesided' forces the return of a one-sided spectrum, while 'twosided' forces two-sided. |
+
+### Returns
+| Variable | Description |
+| --- | --- |
+| Pxx | 1-D array. The values for the power spectrum before scaling (real valued). |
+| freqs | 1-D array. The frequencies corresponding to the elements in `Pxx`. |
+
+### Example Usage
+```python
+from filterbank.filterbank import Filterbank
+import matplotlib.pyplot as plt
+import numpy as np
+from plot import psd
+from filterbank.header import read_header
+
+# Instatiate the filterbank reader and point to the filterbank file
+fb = Filterbank(filename='examples/pspm32.fil')
+
+# read the data in the filterbank file
+_, samples = fb.select_data()
+
+# Convert 2D Array to 1D Array with complex numbers
+samples = samples[0] + (samples[1] * 1j)
+
+# Read the header of the filterbank file
+header = read_header('examples/pspm32.fil')
+
+# Calculate the center frequency with the data in the header
+center_freq = header[b'fch1'] + float(header[b'nchans']) * header[b'foff'] / 2.0
+
+print(center_freq)
+# Get the powerlevels and the frequencies
+PXX, freqs, _ = psd(samples, nfft=1024, sample_rate=80, sides='twosided')
+
+# Calculate the powerlevel dB's
+power_levels =  10 * np.log10(PXX/(80))
+
+# Add the center frequency to the frequencies so they match the actual frequencies
+freqs = freqs + center_freq
+
+# Plot the PSD
+plt.grid(True)
+plt.xlabel('Frequency (MHz)')
+plt.ylabel('Intensity (dB)')
+plt.plot(freqs, power_levels)
+plt.show()
+```

examples/visualize_filterbank_psd.py

@@ -0,0 +1,46 @@
+"""
+    Example of plotting a Power Spectral Density plot, using filterbank data
+"""
+# pylint: disable-all
+import os,sys,inspect
+currentdir = os.path.dirname(os.path.abspath(inspect.getfile(inspect.currentframe())))
+parentdir = os.path.dirname(currentdir)
+sys.path.insert(0,parentdir)
+from filterbank.filterbank import Filterbank
+import matplotlib.pyplot as plt
+import numpy as np
+from plot import psd
+from filterbank.header import read_header
+
+# Instatiate the filterbank reader and point to the filterbank file
+fb = Filterbank(filename='examples/pspm32.fil')
+
+# read the data in the filterbank file
+_, samples = fb.select_data()
+
+# Convert 2D Array to 1D Array with complex numbers
+samples = samples[0] + (samples[1] * 1j)
+
+# Read the header of the filterbank file
+header = read_header('examples/pspm32.fil')
+
+# Calculate the center frequency with the data in the header
+center_freq = header[b'fch1'] + float(header[b'nchans']) * header[b'foff'] / 2.0
+
+print(center_freq)
+# Get the powerlevels and the frequencies
+PXX, freqs = psd(samples, nfft=1024, sample_rate=80,
+                    scale_by_freq=True, sides='twosided')
+
+# Calculate the powerlevel dB's
+power_levels =  10 * np.log10(PXX/(80))
+
+# Add the center frequency to the frequencies so they match the actual frequencies
+freqs = freqs + center_freq
+
+# Plot the PSD
+plt.grid(True)
+plt.xlabel('Frequency (MHz)')
+plt.ylabel('Intensity (dB)')
+plt.plot(freqs, power_levels)
+plt.show()

examples/visualize_psd.py

@@ -26,14 +26,15 @@
 # Close RTLSDR device connection
 SDR.close()
 
+print(SAMPLES[1:100])
 # Number of samples equals the length of samples
 N = SAMPLES.shape[0]
 
 # T equals N/Fs
 T = N/SDR.sample_rate
 
 # Get the powerlevels and the frequencies
-PXX, freqs, _ = psd(SAMPLES, nfft=1024, sample_rate=SDR.sample_rate/1e6, # pylint: disable-msg=C0103
+PXX, freqs = psd(SAMPLES, nfft=1024, sample_rate=SDR.sample_rate/1e6, # pylint: disable-msg=C0103
                     scale_by_freq=True, sides='twosided')
 
 # Calculate the powerlevel dB's

filterbank/__init__.pyc

@@ -7,7 +7,17 @@
 
 def dft_slow(input_data):
     """
-        Compute the discrete Fourier Transform of the one-dimensional array x
+        Compute the discrete Fourier Transform of the one-dimensional array input_data
+
+        Parameters
+        ----------
+        input_data : ndarray
+            Array of length `n` containing the values to be transformed.
+
+        Returns
+        -------
+        ndarray
+            Array of length `n` containing the transformed values.
     """
     input_data = np.asarray(input_data)
     data_length = input_data.shape[0]
@@ -20,6 +30,16 @@ def dft_slow(input_data):
 def fft_vectorized(input_data):
     """
         A vectorized, non-recursive version of the Cooley-Tukey FFT
+
+        Parameters
+        ----------
+        input_data : ndarray
+            Array of length `n` containing the values to be transformed.
+
+        Returns
+        -------
+        ndarray
+            Array of length `n` containing the transformed values.
     """
     input_data = np.asarray(input_data)
     data_length = input_data.shape[0]