Support timestamping in the waterfall intermediate file

The waterfall file has now a constant sized header of 52 bytes,
so that plotting tools can reconstruct properly the spectrum.

The structure of the header is the following:
  - A 32 byte string containing the timestamp in
    ISO-8601 format. This timer has microsecond accuracy.
  - A 4 byte integer containing the sampling rate
  - A 4 byte integer with the FFT size
  - A 4 byte integer containing the number of FFT snapshots for one row
    at the waterfall
  - A 4 byte float with the center frequency of the observation.
  - A 4 byte integer indicating the endianness of the rest of the file. If
    set to 0 the file continues in Big endian. Otherwise, in little endian.
    The change of the endianness is performed to reduce the overhead at the
    station.

 Note that all contents of the header are in Network Byte order! The rest
 of the file is in native byte order, mainly for performance reasons.
 Users can use data of the header to determine if their architecture match
 the architecture of the host generated the waterfall file and act
 accordingly.

 The file continues with information regarding the spectral content of the
 observation.
 Each waterfall line is prepended with a int64_t field indicating the
 absolute time in microseconds with respect to the start of the waterfall
 data (stored in the corresponding header field).
 The spectral content is stored in $FFT$ float values already converted in
 dB scale.
This commit is contained in:
Manolis Surligas 2019-12-16 19:14:33 +02:00
parent 48e421e3f5
commit 8072219a8a
3 changed files with 151 additions and 63 deletions

View File

@ -2,7 +2,7 @@
/*
* gr-satnogs: SatNOGS GNU Radio Out-Of-Tree Module
*
* Copyright (C) 2017, Libre Space Foundation <http://librespacefoundation.org/>
* Copyright (C) 2017,2019 Libre Space Foundation <http://libre.space>
*
* This program is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
@ -31,8 +31,8 @@ namespace satnogs {
* \brief This block computes the waterfall of the incoming signal
* and stores the result to a file.
*
* The file has a special header, so that the satnogs_waterfall Gnuplot
* script to be able to plot it properly.
* The file has a special header, so plotting tools can reconstruct properly
* the spectrum.
*
* \ingroup satnogs
*
@ -45,25 +45,51 @@ public:
* This block computes the waterfall of the incoming signal
* and stores the result to a file.
*
* The file has a special header, so that the satnogs_waterfall Gnuplot
* script to be able to plot it properly.
* The file has a constant sized header of 52 bytes, so that plotting tools can
* reconstruct properly the spectrum.
*
* The structure of the header is the following:
* - A 32 byte string containing the timestamp in
* ISO-8601 format. This timer has microsecond accuracy.
* - A 4 byte integer containing the sampling rate
* - A 4 byte integer with the FFT size
* - A 4 byte integer containing the number of FFT snapshots for one row
* at the waterfall
* - A 4 byte float with the center frequency of the observation.
* - A 4 byte integer indicating the endianness of the rest of the file. If
* set to 0 the file continues in Big endian. Otherwise, in little endian.
* The change of the endianness is performed to reduce the overhead at the
* station.
*
* @note All contents of the header are in Network Byte order! The rest
* of the file is in native byte order, mainly for performance reasons.
* Users can use data of the header to determine if their architecture match
* the architecture of the host generated the waterfall file and act
* accordingly.
*
* The file continues with information regarding the spectral content of the
* observation.
* Each waterfall line is prepended with a int64_t field indicating the
* absolute time in microseconds with respect to the start of the waterfall
* data (stored in the corresponding header field).
* The spectral content is stored in $FFT$ float values already converted in
* dB scale.
*
* @param samp_rate the sampling rate
* @param center_freq the observation center frequency. Used only for
* plotting reasons. For a normalized frequency x-axis set it to 0.
* @param pps pixels per second
* @param rps rows per second
* @param fft_size FFT size
* @param filename the name of the output file
* @param mode the mode that the waterfall.
* - 0: Simple decimation
* - 1: Max hold
* - 2: Mean energy
*
* @return shared pointer to the object
*/
static sptr
make(double samp_rate, double center_freq,
double pps, size_t fft_size,
make(float samp_rate, float center_freq,
float rps, size_t fft_size,
const std::string &filename, int mode = 0);
};

View File

@ -2,7 +2,7 @@
/*
* gr-satnogs: SatNOGS GNU Radio Out-Of-Tree Module
*
* Copyright (C) 2017, Libre Space Foundation <http://librespacefoundation.org/>
* Copyright (C) 2017,2019 Libre Space Foundation <http://libre.space>
*
* This program is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
@ -25,56 +25,55 @@
#include <gnuradio/io_signature.h>
#include "waterfall_sink_impl.h"
#include <satnogs/log.h>
#include <satnogs/utils.h>
#include <satnogs/date.h>
namespace gr {
namespace satnogs {
waterfall_sink::sptr
waterfall_sink::make(double samp_rate, double center_freq,
double fps, size_t fft_size,
const std::string &filename, int mode)
waterfall_sink::make(float samp_rate, float center_freq, float rps,
size_t fft_size, const std::string &filename, int mode)
{
return gnuradio::get_initial_sptr(
new waterfall_sink_impl(samp_rate, center_freq,
fps, fft_size, filename, mode));
new waterfall_sink_impl(samp_rate, center_freq, rps, fft_size, filename,
mode));
}
/*
* The private constructor
*/
waterfall_sink_impl::waterfall_sink_impl(double samp_rate,
double center_freq,
double pps,
size_t fft_size,
waterfall_sink_impl::waterfall_sink_impl(float samp_rate, float center_freq,
float rps, size_t fft_size,
const std::string &filename,
int mode) :
gr::sync_block("waterfall_sink",
gr::io_signature::make(1, 1, sizeof(gr_complex)),
gr::io_signature::make(0, 0, 0)),
d_samp_rate(samp_rate),
d_pps(pps),
d_center_freq(center_freq),
d_fft_size(fft_size),
d_mode((wf_mode_t)mode),
d_refresh((d_samp_rate / fft_size) / pps),
d_mode((wf_mode_t) mode),
d_refresh((d_samp_rate / fft_size) / rps),
d_fft_cnt(0),
d_fft_shift((size_t)(ceil(fft_size / 2.0))),
d_samples_cnt(0),
d_fft(fft_size)
{
float r = 0.0;
const int alignment_multiple = volk_get_alignment()
/ (fft_size * sizeof(gr_complex));
set_alignment(std::max(1, alignment_multiple));
set_output_multiple(fft_size);
d_shift_buffer = (gr_complex *) volk_malloc(
fft_size * sizeof(gr_complex), volk_get_alignment());
d_shift_buffer = (gr_complex *) volk_malloc(fft_size * sizeof(gr_complex),
volk_get_alignment());
if (!d_shift_buffer) {
LOG_ERROR("Could not allocate aligned memory");
throw std::runtime_error("Could not allocate aligned memory");
}
d_hold_buffer = (float *)volk_malloc(fft_size * sizeof(gr_complex),
volk_get_alignment());
d_hold_buffer = (float *) volk_malloc(fft_size * sizeof(gr_complex),
volk_get_alignment());
if (!d_hold_buffer) {
LOG_ERROR("Could not allocate aligned memory");
throw std::runtime_error("Could not allocate aligned memory");
@ -89,16 +88,23 @@ waterfall_sink_impl::waterfall_sink_impl(double samp_rate,
}
d_fos.open(filename, std::ios::binary | std::ios::trunc);
/* Append header for proper plotting */
r = fft_size;
d_fos.write((char *)&r, sizeof(float));
for (size_t i = 0; i < fft_size; i++) {
r = (samp_rate / fft_size * i) - samp_rate / 2.0 + center_freq;
d_fos.write((char *)&r, sizeof(float));
if (d_fos.fail()) {
throw std::runtime_error("Could not create file for writing");
}
}
bool
waterfall_sink_impl::start()
{
/*
* Append header for proper plotting. We do it on the start() to reduce
* as much as possible the delay between the start of the observation tagging
* and the fist invocation of the work() method.
*/
apply_header();
return true;
}
/*
* Our virtual destructor.
*/
@ -133,7 +139,6 @@ waterfall_sink_impl::work(int noutput_items,
throw std::runtime_error("Wrong waterfall mode");
return -1;
}
return n_fft * d_fft_size;
}
@ -141,11 +146,10 @@ void
waterfall_sink_impl::compute_decimation(const gr_complex *in, size_t n_fft)
{
size_t i;
float t;
gr_complex *fft_in;
for (i = 0; i < n_fft; i++) {
d_fft_cnt++;
if (d_fft_cnt > d_refresh) {
if (d_fft_cnt == d_refresh) {
fft_in = d_fft.get_inbuf();
memcpy(fft_in, in + i * d_fft_size, d_fft_size * sizeof(gr_complex));
d_fft.execute();
@ -161,8 +165,7 @@ waterfall_sink_impl::compute_decimation(const gr_complex *in, size_t n_fft)
(float) d_fft_size, 1.0,
d_fft_size);
/* Write the result to the file */
t = (float)(d_samples_cnt / d_samp_rate);
d_fos.write((char *) &t, sizeof(float));
write_timestamp();
d_fos.write((char *) d_hold_buffer, d_fft_size * sizeof(float));
d_fft_cnt = 0;
}
@ -175,7 +178,6 @@ waterfall_sink_impl::compute_max_hold(const gr_complex *in, size_t n_fft)
{
size_t i;
size_t j;
float t;
gr_complex *fft_in;
for (i = 0; i < n_fft; i++) {
fft_in = d_fft.get_inbuf();
@ -184,29 +186,27 @@ waterfall_sink_impl::compute_max_hold(const gr_complex *in, size_t n_fft)
/* Perform FFT shift */
memcpy(d_shift_buffer, &d_fft.get_outbuf()[d_fft_shift],
sizeof(gr_complex) * (d_fft_size - d_fft_shift));
memcpy(&d_shift_buffer[d_fft_size - d_fft_shift],
&d_fft.get_outbuf()[0], sizeof(gr_complex) * d_fft_shift);
memcpy(&d_shift_buffer[d_fft_size - d_fft_shift], &d_fft.get_outbuf()[0],
sizeof(gr_complex) * d_fft_shift);
/* Normalization factor */
volk_32fc_s32fc_multiply_32fc(d_shift_buffer, d_shift_buffer,
1.0 / d_fft_size, d_fft_size);
/* Compute the mag^2 */
volk_32fc_magnitude_squared_32f(d_tmp_buffer, d_shift_buffer,
d_fft_size);
volk_32fc_magnitude_squared_32f(d_tmp_buffer, d_shift_buffer, d_fft_size);
/* Max hold */
volk_32f_x2_max_32f(d_hold_buffer, d_hold_buffer, d_tmp_buffer,
d_fft_size);
d_fft_cnt++;
if (d_fft_cnt > d_refresh) {
if (d_fft_cnt == d_refresh) {
/* Compute the energy in dB */
for (j = 0; j < d_fft_size; j++) {
d_hold_buffer[j] = 10.0 * log10f(d_hold_buffer[j] + 1.0e-20);
}
/* Write the result to the file */
t = (float)(d_samples_cnt / d_samp_rate);
d_fos.write((char *) &t, sizeof(float));
write_timestamp();
d_fos.write((char *) d_hold_buffer, d_fft_size * sizeof(float));
/* Reset */
@ -217,12 +217,51 @@ waterfall_sink_impl::compute_max_hold(const gr_complex *in, size_t n_fft)
}
}
void
waterfall_sink_impl::apply_header()
{
header_t h;
memset(h.start_time, 0, 32);
std::chrono::system_clock::time_point tp = std::chrono::system_clock::now();
d_start = tp;
std::string s = date::format("%FT%TZ",
date::floor<std::chrono::microseconds> (tp));
std::strncpy(h.start_time, s.c_str(), 32);
/* Before writing to the file convert all values to Network Byte Order */
h.fft_size = htonl(d_fft_size);
h.samp_rate = htonl(d_samp_rate);
h.nfft_per_row = htonl(d_refresh);
uint32_t tmp = htonl(*((uint32_t *) &d_center_freq));
memcpy(&h.center_freq, &tmp, sizeof(uint32_t));
h.center_freq = *((float *) &tmp);
h.endianness = !(1 == htonl(1));
/*
* Write the header. Make a dummy serialization to avoid padding and
* alignment issues
*/
d_fos.write(h.start_time, 32);
d_fos.write((char *)&h.fft_size, sizeof(uint32_t));
d_fos.write((char *)&h.samp_rate, sizeof(uint32_t));
d_fos.write((char *)&h.nfft_per_row, sizeof(uint32_t));
d_fos.write((char *)&h.center_freq, sizeof(float));
d_fos.write((char *)&h.endianness, sizeof(uint32_t));
}
void
waterfall_sink_impl::write_timestamp()
{
std::chrono::system_clock::time_point tp = std::chrono::system_clock::now();
int64_t x = std::chrono::duration_cast<std::chrono::microseconds> (
tp - d_start).count();
d_fos.write((char *)&x, sizeof(int64_t));
}
void
waterfall_sink_impl::compute_mean(const gr_complex *in, size_t n_fft)
{
size_t i;
size_t j;
float t;
gr_complex *fft_in;
for (i = 0; i < n_fft; i++) {
fft_in = d_fft.get_inbuf();
@ -231,25 +270,24 @@ waterfall_sink_impl::compute_mean(const gr_complex *in, size_t n_fft)
/* Perform FFT shift */
memcpy(d_shift_buffer, &d_fft.get_outbuf()[d_fft_shift],
sizeof(gr_complex) * (d_fft_size - d_fft_shift));
memcpy(&d_shift_buffer[d_fft_size - d_fft_shift],
&d_fft.get_outbuf()[0], sizeof(gr_complex) * d_fft_shift);
memcpy(&d_shift_buffer[d_fft_size - d_fft_shift], &d_fft.get_outbuf()[0],
sizeof(gr_complex) * d_fft_shift);
/* Accumulate the complex numbers */
volk_32f_x2_add_32f(d_hold_buffer, d_hold_buffer,
(float *)d_shift_buffer, 2 * d_fft_size);
volk_32f_x2_add_32f(d_hold_buffer, d_hold_buffer, (float *) d_shift_buffer,
2 * d_fft_size);
d_fft_cnt++;
if (d_fft_cnt > d_refresh) {
if (d_fft_cnt == d_refresh) {
/*
* Compute the energy in dB performing the proper normalization
* before any dB calculation, emulating the mean
*/
volk_32fc_s32f_x2_power_spectral_density_32f(
d_hold_buffer, (gr_complex *)d_hold_buffer,
d_hold_buffer, (gr_complex *) d_hold_buffer,
(float) d_fft_cnt * d_fft_size, 1.0, d_fft_size);
/* Write the result to the file */
t = (float)(d_samples_cnt / d_samp_rate);
d_fos.write((char *) &t, sizeof(float));
write_timestamp();
d_fos.write((char *) d_hold_buffer, d_fft_size * sizeof(float));
/* Reset */

View File

@ -2,7 +2,7 @@
/*
* gr-satnogs: SatNOGS GNU Radio Out-Of-Tree Module
*
* Copyright (C) 2017, Libre Space Foundation <http://librespacefoundation.org/>
* Copyright (C) 2017,2019 Libre Space Foundation <http://libre.space>
*
* This program is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
@ -26,12 +26,28 @@
#include <gnuradio/fft/fft.h>
#include <iostream>
#include <fstream>
#include <chrono>
namespace gr {
namespace satnogs {
class waterfall_sink_impl : public waterfall_sink {
private:
/**
* Waterfall header data.
* This structure is only for readability purposes and make more clear to
* possible users the structure of the header.
*/
typedef struct {
char start_time[32]; /**< String with the start of the waterfall in ISO-8601 format */
uint32_t samp_rate; /**< The sampling rate of the waterfall */
uint32_t fft_size; /**< The FFT size of the flowgraph */
uint32_t nfft_per_row; /**< The number of FFTs performed to plot one row */
float center_freq; /**< The center frequency of the observation. Just for viasualization purposes */
uint32_t endianness; /**< The endianness of the rest of the file. Should be 0 for big endian */
} header_t;
/**
* The different types of operation of the waterfall
*/
@ -41,8 +57,8 @@ private:
WATERFALL_MODE_MEAN = 2 //!< WATERFALL_MODE_MEAN compute the mean energy of all the FFT snapshots between two consecutive pixel rows
} wf_mode_t;
const double d_samp_rate;
double d_pps;
const float d_samp_rate;
const float d_center_freq;
const size_t d_fft_size;
wf_mode_t d_mode;
size_t d_refresh;
@ -54,13 +70,21 @@ private:
float *d_hold_buffer;
float *d_tmp_buffer;
std::ofstream d_fos;
std::chrono::system_clock::time_point d_start;
void
apply_header();
void
write_timestamp();
public:
waterfall_sink_impl(double samp_rate, double center_freq,
double pps, size_t fft_size,
const std::string &filename, int mode);
waterfall_sink_impl(float samp_rate, float center_freq, float rps,
size_t fft_size, const std::string &filename, int mode);
~waterfall_sink_impl();
bool
start();
int
work(int noutput_items, gr_vector_const_void_star &input_items,