SSB_HighSpeed_Modem/hsmodem/fft.cpp

188 lines
5.0 KiB
C++
Raw Normal View History

2020-11-05 13:11:57 -05:00
/*
* High Speed modem to transfer data in a 2,7kHz SSB channel
* =========================================================
* Author: DJ0ABR
*
* (c) DJ0ABR
* www.dj0abr.de
*
* 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
* the Free Software Foundation; either version 2 of the License, or
* (at your option) any later version.
*
* This program is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program; if not, write to the Free Software
* Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
*
*/
#include "hsmodem.h"
#ifdef _WIN32_
#include "fftw_lib/fftw3.h"
#endif
#ifdef _LINUX_
#include <fftw3.h>
#endif
#define FFT_AUDIOSAMPLERATE 8000
double *din = NULL; // input data for fft
fftw_complex *cpout = NULL; // ouput data from fft
fftw_plan plan = NULL;
#define fft_rate (FFT_AUDIOSAMPLERATE / 10) // resolution: 10 Hz
int fftidx = 0;
int fftcnt = fft_rate/2+1; // number of output values
uint16_t fftout[FFT_AUDIOSAMPLERATE / 10/2+1];
2020-11-24 10:07:52 -05:00
float f_fftout[FFT_AUDIOSAMPLERATE / 10 / 2 + 1];
2020-11-05 13:11:57 -05:00
int downsamp = 0;
int downphase = 0;
2020-11-24 10:07:52 -05:00
int rxlevel_deteced = 0;
int rx_in_sync = 0;
2020-11-05 13:11:57 -05:00
uint16_t *make_waterfall(float fre, int *retlen)
{
// Downsampling:
// needed 8000 bit/s
// caprate 48k: downsample by 6
// caprate 44,1k: downsample by 5,5
2020-12-11 16:53:26 -05:00
if (physcaprate == 48000)
2020-11-05 13:11:57 -05:00
if (++downsamp < 6) return NULL;
2020-11-24 10:07:52 -05:00
2020-12-11 16:53:26 -05:00
// TODO: the following simple resamp results in double peeks in fft
if (physcaprate == 44100)
2020-11-05 13:11:57 -05:00
{
if (downphase <= 1100)
{
if (++downsamp < 5) return NULL;
}
else
{
if (++downsamp < 6) return NULL;
}
if(++downphase >= 2000) downphase = 0;
}
downsamp = 0;
int fftrdy = 0;
// fre are the float samples
// fill into the fft input buffer
din[fftidx++] = fre;
if(fftidx == fft_rate)
{
fftidx = 0;
2020-11-24 10:07:52 -05:00
2020-11-05 13:11:57 -05:00
// the fft buffer is full, execute the FFT
fftw_execute(plan);
for (int j = 0; j < fftcnt; j++)
{
// calculate absolute value (magnitute without phase)
float fre = (float)cpout[j][0];
float fim = (float)cpout[j][1];
float mag = sqrt((fre * fre) + (fim * fim));
fftout[j] = (uint16_t)mag;
2020-11-24 10:07:52 -05:00
f_fftout[j] = mag;
2020-11-05 13:11:57 -05:00
fftrdy = 1;
}
2020-11-24 10:07:52 -05:00
2020-12-11 16:53:26 -05:00
if (rx_in_sync == 0)
2020-11-24 10:07:52 -05:00
{
2020-12-11 16:53:26 -05:00
// signal detection
// measure level at band edges
float edgelevel = 0;
for (int e = 0; e < 10; e++)
edgelevel += f_fftout[e];
edgelevel /= 10;
for (int e = 300; e < 320; e++)
edgelevel += f_fftout[e];
edgelevel /= 20;
// measure level at mid band
float midlevel = 0;
for (int e = 100; e < 300; e++)
midlevel += f_fftout[e];
midlevel /= 200;
//calc difference in %
int idiff = (int)((edgelevel * 100) / midlevel);
//printf("diff:%d %% edge:%10.6f midband:%10.6f\n", idiff,edgelevel, midlevel);
// idiff SDR Console:
// no signal ... > 100
// signal < 20
static int checks = 0;
static int lastsig = 0;
int sig = 0;
// check if signal detected or not
if (idiff > 100) sig = 0;
if (idiff < 20) sig = 1;
rxlevel_deteced = sig;
// check if changed since last check
if (sig != lastsig)
2020-11-24 10:07:52 -05:00
{
2020-12-11 16:53:26 -05:00
lastsig = sig;
checks = 0;
}
else
{
if (checks <= 3) checks++;
if (checks == 3)
{
if (sig == 1)
{
printf("===>>> level detected, reset modem\n");
trigger_resetmodem = 1;
}
}
2020-11-24 10:07:52 -05:00
}
}
2020-11-05 13:11:57 -05:00
}
if(fftrdy == 1)
{
*retlen = fftcnt;
return fftout;
}
return NULL;
}
void init_fft()
{
2020-11-21 13:54:47 -05:00
/*
2020-11-05 13:11:57 -05:00
char fn[300];
2020-11-21 13:54:47 -05:00
* storing to a file in the working directory may be a problem under Windows, so we do not use wisdom files
sprintf(fn, "capture_fft_%d", fft_rate); // wisdom file for each capture rate
2020-11-05 13:11:57 -05:00
2020-11-21 13:54:47 -05:00
fftw_import_wisdom_from_filename(fn);*/
2020-11-05 13:11:57 -05:00
din = (double *)fftw_malloc(sizeof(double) * fft_rate);
cpout = (fftw_complex *)fftw_malloc(sizeof(fftw_complex) * fft_rate);
plan = fftw_plan_dft_r2c_1d(fft_rate, din, cpout, FFTW_MEASURE);
2020-11-21 13:54:47 -05:00
//fftw_export_wisdom_to_filename(fn);
2020-11-05 13:11:57 -05:00
}
void exit_fft()
{
if(plan) fftw_destroy_plan(plan);
if(din) fftw_free(din);
if(cpout) fftw_free(cpout);
}