diff options
| -rw-r--r-- | doc/example.ini | 13 | ||||
| -rwxr-xr-x | doc/time-freq-plot.py | 77 | ||||
| -rw-r--r-- | src/DabMod.cpp | 14 | ||||
| -rw-r--r-- | src/Resampler.cpp | 4 | 
4 files changed, 101 insertions, 7 deletions
| diff --git a/doc/example.ini b/doc/example.ini index 0659813..f8cec36 100644 --- a/doc/example.ini +++ b/doc/example.ini @@ -141,13 +141,16 @@ output=uhd  [fileoutput]  ; Two output formats are supported: In the default mode, -; the file output writes I/Q float values (i.e. complex float) -; to the file. The I and Q samples can take values up to -; 100000 in absolute magnitude with gainmode FIX. -; With gainmode VAR, they should never exceed 50000. -; With gainmode MAX, thet are limited to 32767. +; the file output writes I/Q float values (i.e. complex +; float) to the file. The I and Q samples can take values up +; to 810000 in absolute magnitude with gainmode FIX. With +; gainmode VAR and FIX, they should never exceed 50000.  ;format=complexf  ; +; When the format is set to complexf_normalised the maximal +; amplitude written to a file, is the digital gain +-10%. +;format=complexf_normalised +;  ; When the format is set to s8, the output writes I/Q 8-bit  ; signed integers, where the magnitude is multiplied by 128/50000  ; effectively mapping the gainmode VAR range of -50000 -- 50000 diff --git a/doc/time-freq-plot.py b/doc/time-freq-plot.py new file mode 100755 index 0000000..6ece564 --- /dev/null +++ b/doc/time-freq-plot.py @@ -0,0 +1,77 @@ +#!/usr/bin/python +# -*- coding: utf-8 -*- +# +# Print scope and spectrum from ODR-DabMod I/Q file +# +# The MIT License (MIT) +# +# Copyright (c) 2017 Matthias P. Braendli +# +# Permission is hereby granted, free of charge, to any person obtaining a copy +# of this software and associated documentation files (the "Software"), to deal +# in the Software without restriction, including without limitation the rights +# to use, copy, modify, merge, publish, distribute, sublicense, and/or sell +# copies of the Software, and to permit persons to whom the Software is +# furnished to do so, subject to the following conditions: +# +#    The above copyright notice and this permission notice shall be included in +#    all copies or substantial portions of the Software. +# +#    THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR +#    IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, +#    FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE +#    AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER +#    LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, +#    OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN +#    THE SOFTWARE. +import sys +import matplotlib.pyplot as plt +import numpy as np + +rate=2048000 +# T = 1/2048000 s +# NULL symbol is 2656 T (about 1.3ms) long. +T_NULL = 2656 +# Full transmission frame in TM1 is 96ms = 196608 T. +T_TF = 196608 + +num_skip_samples = 8 * T_TF +num_analyse_samples = 2 * T_TF + +if len(sys.argv) < 2: +    print("Specify .iq file name") +    print("Expected format: complex float I/Q, 2048000 Sps") +    print("The input file must contain at least 10 transmission frames,") +    print("i.e. {} samples = {} seconds".format(T_TF * 10, T_TF * 10.0 / rate)) +    sys.exit(1) + +fd = open(sys.argv[1], 'rb') + +# The IQ files potentially have zero samples in the beginning, we need +# to skip a few transmission frames + + +source_data = np.fromfile(file=fd, dtype=np.complex64, count=num_skip_samples + num_analyse_samples) + +print("Read in {} samples".format(len(source_data))) + +source_data = source_data[num_skip_samples:] +source_data_time = np.linspace(0, num_analyse_samples/rate, len(source_data)) + +print("Signal power: {} of {} samples".format(np.sum(np.abs(source_data**2)), len(source_data))) + +fft_size = 4096 + +plt.figure(figsize=(10,8)) +plt.subplot(211) +plt.title("Real part of signal") +plt.plot(source_data_time, np.real(source_data)) + +signal_spectrum = np.abs(np.fft.fftshift(np.fft.fft(source_data[T_NULL:], fft_size))) +freqs = np.fft.fftshift(np.fft.fftfreq(fft_size, d=1./rate)) + +plt.subplot(212) +plt.title("Spectrum of {} samples after the NULL symbol".format(fft_size)) +plt.semilogy(freqs, signal_spectrum) + +plt.show() diff --git a/src/DabMod.cpp b/src/DabMod.cpp index 301b078..4e4cdab 100644 --- a/src/DabMod.cpp +++ b/src/DabMod.cpp @@ -81,6 +81,11 @@   */  static const float normalise_factor = 50000.0f; +//Empirical normalisation factors used to normalise the samples to amplitude 1. +static const float normalise_factor_file_fix = 81000.0f; +static const float normalise_factor_file_var = 46000.0f; +static const float normalise_factor_file_max = 46000.0f; +  typedef std::complex<float> complexf;  using namespace std; @@ -183,6 +188,15 @@ static shared_ptr<ModOutput> prepare_output(          if (s.fileOutputFormat == "complexf") {              output = make_shared<OutputFile>(s.outputName);          } +        if (s.fileOutputFormat == "complexf_normalised") { +            if (s.gainMode == GainMode::GAIN_FIX) +                s.normalise = 1.0 / normalise_factor_file_fix; +            else if (s.gainMode == GainMode::GAIN_MAX) +                s.normalise = 1.0 / normalise_factor_file_max; +            else if (s.gainMode == GainMode::GAIN_VAR) +                s.normalise = 1.0 / normalise_factor_file_var; +            output = make_shared<OutputFile>(s.outputName); +        }          else if (s.fileOutputFormat == "s8") {              // We must normalise the samples to the interval [-127.0; 127.0]              s.normalise = 127.0f / normalise_factor; diff --git a/src/Resampler.cpp b/src/Resampler.cpp index ee2b865..8786e91 100644 --- a/src/Resampler.cpp +++ b/src/Resampler.cpp @@ -76,9 +76,9 @@ Resampler::Resampler(size_t inputRate, size_t outputRate, size_t resolution) :      PDEBUG(" FFT size in: %zu, FFT size out: %zu\n", myFftSizeIn, myFftSizeOut);      if (myFftSizeIn > myFftSizeOut) { -        myFactor = 1.0f / myFftSizeIn; +        myFactor = 1.0f / myFftSizeIn * outputRate / inputRate;      } else { -        myFactor = 1.0f / myFftSizeOut; +        myFactor = 1.0f / myFftSizeOut * outputRate / inputRate;      }      myWindow = (float*)memalign(16, myFftSizeIn * sizeof(float)); | 
