/* * osmo-fl2k, turns FL2000-based USB 3.0 to VGA adapters into * low cost DACs * * fl2k-iq * Copyright (C) 2020 by Felix Erckenbrecht * * based on fl2k-fm code: * Copyright (C) 2016-2018 by Steve Markgraf * * based on FM modulator code from VGASIG: * Copyright (C) 2009 by Bartek Kania * * SPDX-License-Identifier: GPL-2.0+ * * 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, see . */ #include #include #include #include #include #ifndef _WIN32 #include #include #include #else #include #include #include #include "getopt/getopt.h" #endif #include #include #include #include "osmo-fl2k.h" enum inputType_E { INP_REAL, INP_COMPLEX }; #define BUFFER_SAMPLES_SHIFT 16 #define BUFFER_SAMPLES (1 << BUFFER_SAMPLES_SHIFT) #define BUFFER_SAMPLES_MASK ((1 << BUFFER_SAMPLES_SHIFT)-1) #define BASEBAND_BUF_SIZE 2048 fl2k_dev_t *dev = NULL; volatile int do_exit = 0; pthread_t iq_thread; pthread_mutex_t cb_mutex; pthread_mutex_t iq_mutex; pthread_cond_t cb_cond; pthread_cond_t iq_cond; FILE *file; int8_t *itxbuf = NULL; int8_t *qtxbuf = NULL; int8_t *iambuf = NULL; int8_t *qambuf = NULL; int8_t *buf1 = NULL; int8_t *buf2 = NULL; uint32_t samp_rate = 96000000; int base_freq = 1440000; int rf_to_baseband_sample_ratio; int input_freq = 48000; complex float *ampbuf; complex float *slopebuf; long int * pdbuf; long int * pdslopebuf; int writepos, readpos; int swap_iq = 0; int ignore_eof = 0; void usage(void) { fprintf(stderr, "fl2k_ampliphase, a special modulator for FL2K VGA dongles\n\n" "Usage:" "\t[-d device index (default: 0)]\n" "\t[-c center frequency (default: 1440 kHz)]\n" "\t[-i input baseband sample rate (default: 48000 Hz)]\n" "\t[-s samplerate in Hz (default: 96 MS/s)]\n" "\t[-t type of input: real/complex (default: real)\n" "\t input requirements: real - single channel (mono)\n" "\t complex - dual channel (stereo)\n" "\t[-w swap I & Q (invert spectrum)\n" "\t[-e ignore EOF\n" "\tfilename (use '-' to read from stdin)\n\n" ); exit(1); } #ifdef _WIN32 BOOL WINAPI sighandler(int signum) { if (CTRL_C_EVENT == signum) { fprintf(stderr, "Signal caught, exiting!\n"); fl2k_stop_tx(dev); do_exit = 1; pthread_cond_signal(&iq_cond); return TRUE; } return FALSE; } #else static void sighandler(int signum) { fprintf(stderr, "Signal caught, exiting!\n"); fl2k_stop_tx(dev); do_exit = 1; pthread_cond_signal(&iq_cond); } #endif /* DDS Functions */ #ifndef M_PI # define M_PI 3.14159265358979323846 /* pi */ # define M_PI_2 1.57079632679489661923 /* pi/2 */ # define M_PI_4 0.78539816339744830962 /* pi/4 */ # define M_1_PI 0.31830988618379067154 /* 1/pi */ # define M_2_PI 0.63661977236758134308 /* 2/pi */ #endif #define DDS_2PI (M_PI * 2) /* 2 * Pi */ #define DDS_3PI2 (M_PI_2 * 3) /* 3/2 * pi */ #define TRIG_TABLE_ORDER 8 #define TRIG_TABLE_SHIFT (32 - TRIG_TABLE_ORDER) #define TRIG_TABLE_LEN (1 << TRIG_TABLE_ORDER) //#define ANG_INCR (0xffffffff / DDS_2PI) #define ANG_INCR ((float)(0x100000000)) / DDS_2PI enum waveform_E { WF_SINE, WF_RECT }; struct trigonometric_table_S { int initialized; int16_t quadrature[TRIG_TABLE_LEN]; int16_t inphase[TRIG_TABLE_LEN]; }; static struct trigonometric_table_S trig_table = { .initialized = 0 }; typedef struct { float sample_freq; float freq; /* instantaneous phase */ unsigned long int phase; /* phase increment */ unsigned long int phase_step; /* for phase modulation */ long int phase_delta; long int phase_slope; /* for amplitude modulation */ complex float amplitude; complex float ampslope; } dds_t; static inline void dds_set_freq(dds_t *dds, float freq) { dds->freq = freq; dds->phase_step = (freq / dds->sample_freq) * 2 * M_PI * ANG_INCR; } static inline void dds_set_amp(dds_t *dds, complex float amplitude, complex float ampslope) { dds->amplitude = amplitude; dds->ampslope = ampslope; } static inline void dds_set_phase(dds_t *dds, long int phase_delta, long int phase_slope) { dds->phase_delta = phase_delta; dds->phase_slope = phase_slope; } dds_t dds_init(float sample_freq, float freq, float phase, float amp, enum waveform_E waveform ) { dds_t dds; int i; dds.sample_freq = sample_freq; dds.phase = phase * ANG_INCR; dds_set_freq(&dds, freq); dds_set_amp(&dds, amp, 0); /* Initialize quadrature table, prescaled for 16 bit signed integer */ if (!trig_table.initialized) { float incr = 1.0 / (float)TRIG_TABLE_LEN; for (i = 0; i < TRIG_TABLE_LEN; i++){ if(waveform == WF_SINE){ trig_table.quadrature[i]= sin(incr * i * DDS_2PI) * 32767; trig_table.inphase[i] = cos(incr * i * DDS_2PI) * 32767; } else{ /* rectangular / square output */ trig_table.quadrature[i]= sin(incr * i * DDS_2PI) >= 0 ? 32767 : -32767; trig_table.inphase[i] = cos(incr * i * DDS_2PI) >= 0 ? 32767 : -32767; } } trig_table.initialized = 1; } return dds; } static inline int8_t dds_real(dds_t *dds) { int tmp; int32_t amp_i, amp_q; int8_t amp8; // advance dds generator tmp = dds->phase >> TRIG_TABLE_SHIFT; dds->phase += dds->phase_step; dds->phase &= 0xffffffff; //amp = 255; amp_i = creal(dds->amplitude) * 23170.0; // 0..15, * 1/SQRT(2) amp_q = cimag(dds->amplitude) * 23170.0; amp_i = amp_i * trig_table.quadrature[tmp]; // 0..31, * 1/SQRT(2) amp_q = amp_q * trig_table.inphase[tmp]; // 0..31, * 1/SQRT(2) amp8 = (int8_t) ((amp_i + amp_q) >> 24); // 0..31 >> 24 => 0..8 dds->amplitude += dds->ampslope; return amp8; } static inline void dds_real_buf(dds_t *dds, int8_t *buf, int count) { int i; for (i = 0; i < count; i++) buf[i] = dds_real(dds); } static inline void dds_complex(dds_t *dds, int8_t * i, int8_t * q) { int pi_i, pi_q; int32_t amp_i, amp_q; // get current carrier phase, add phase mod, calculate table index pi_i = (dds->phase - dds->phase_delta) >> TRIG_TABLE_SHIFT; pi_q = (dds->phase + dds->phase_delta) >> TRIG_TABLE_SHIFT; // advance dds generator dds->phase += dds->phase_step; // add some extra phase modulation dds->phase += dds->phase_delta; dds->phase &= 0xffffffff; //amp = 255; amp_i = (int32_t) (creal(dds->amplitude) * 32767.0); // 0..15 amp_q = (int32_t) (cimag(dds->amplitude) * 32767.0); amp_i = amp_i * trig_table.inphase[pi_i]; // 0..31 amp_q = amp_q * trig_table.quadrature[pi_q]; // 0..31 *i = (int8_t) (amp_i >> 24); // 0..31 >> 24 => 0..8 *q = (int8_t) (amp_q >> 24); // 0..31 >> 24 => 0..8 dds->amplitude += dds->ampslope; dds->phase_delta += dds->phase_slope; return; } static inline void dds_complex_buf(dds_t *dds, int8_t *ibuf, int8_t *qbuf, int count) { int i; for (i = 0; i < count; i++){ dds_complex(dds, &ibuf[i], &qbuf[i]); } } /* Signal generation and some helpers */ /* Generate the radio signal using the pre-calculated amplitude information * in the amp buffer */ static void *iq_worker(void *arg) { register float freq; register float tmp; dds_t base_signal; int8_t *tmp_ptr; uint32_t len = 0; uint32_t readlen, remaining; int buf_prefilled = 0; /* Prepare the oscillators */ base_signal = dds_init(samp_rate, base_freq, 0, 1, WF_RECT); while (!do_exit) { // dds_set_amp(&base_signal, ampbuf[readpos], slopebuf[readpos]); /* phase modulate the oscillator */ dds_set_phase(&base_signal, pdbuf[readpos], pdslopebuf[readpos]); readpos++; readpos &= BUFFER_SAMPLES_MASK; /* check if we reach the end of the buffer */ if ((len + rf_to_baseband_sample_ratio) > FL2K_BUF_LEN) { readlen = FL2K_BUF_LEN - len; remaining = rf_to_baseband_sample_ratio - readlen; dds_complex_buf(&base_signal, &iambuf[len], &qambuf[len],readlen); if (buf_prefilled) { /* swap buffers */ tmp_ptr = iambuf; iambuf = itxbuf; itxbuf = tmp_ptr; tmp_ptr = qambuf; qambuf = qtxbuf; qtxbuf = tmp_ptr; pthread_mutex_lock(&cb_mutex); pthread_cond_wait(&cb_cond, &cb_mutex); pthread_mutex_unlock(&cb_mutex); } dds_complex_buf(&base_signal, iambuf, qambuf, remaining); len = remaining; buf_prefilled = 1; } else { dds_complex_buf(&base_signal, &iambuf[len], &qambuf[len], rf_to_baseband_sample_ratio); len += rf_to_baseband_sample_ratio; } pthread_mutex_lock(&iq_mutex); pthread_cond_signal(&iq_cond); pthread_mutex_unlock(&iq_mutex); } pthread_exit(NULL); } static inline int writelen(int maxlen) { int rp = readpos; int len; int r; if (rp < writepos) rp += BUFFER_SAMPLES; len = rp - writepos; r = len > maxlen ? maxlen : len; return r; } static inline float complex modulate_sample_ampliphase(const int lastwritepos, const float lastamp, const float sample, float modulationIndex) { float amp; float slope; /* Calculate modulator amplitudes at this point to lessen * the calculations needed in the signal generator */ amp = sample; /* What we do here is calculate a linear "slope" from the previous sample to this one. This is then used by the modulator to gently increase/decrease the phase with each sample without the need to recalculate the dds parameters. In fact this gives us a very efficient and pretty good interpolation filter. */ slope = amp - lastamp; slope = slope * 1.0/ (float) rf_to_baseband_sample_ratio; pdbuf[writepos] = (long int) lastamp * modulationIndex * ANG_INCR; pdslopebuf[writepos] = (long int) slope * modulationIndex * ANG_INCR; return amp; } void ampliphase_modulator(enum inputType_E inputType) { unsigned int i; size_t len; float freq; float complex lastamp = 0; int16_t baseband_buf_real[BASEBAND_BUF_SIZE]; int16_t baseband_buf_cplx[BASEBAND_BUF_SIZE][2]; uint32_t lastwritepos = writepos; float complex sample; while (!do_exit) { int swap = swap_iq; len = writelen(BASEBAND_BUF_SIZE); if (len > 1) { if(inputType == INP_REAL){ len = fread(baseband_buf_real, 1, len, file); for(i = 0 ; i < len; i++){ /* input is -1.0 .. +1.0 (-32768 .. 32767) * transform to 0.0 .. +1.0 (0 .. 32767) * put into I part of BB */ baseband_buf_cplx[i][0] = baseband_buf_real[i] / 2 + INT16_MAX/2; /* Q part of BB is zero for AM */ baseband_buf_cplx[i][1] = 0; } } else{ len = fread(baseband_buf_cplx, 2, len, file); } if (len == 0){ if(ferror(file)){ do_exit = 1; } if(!ignore_eof && feof(file)){ do_exit = 1; } } for (i = 0; i < len; i++) { sample = (float) baseband_buf_cplx[i][0+swap] / 32768.0 + I * (float) baseband_buf_cplx[i][1-swap] / 32768.0; /* Modulate and buffer the sample */ lastamp = modulate_sample_ampliphase(lastwritepos, lastamp, sample); lastwritepos = writepos++; writepos %= BUFFER_SAMPLES; } } else { pthread_mutex_lock(&iq_mutex); pthread_cond_wait(&iq_cond, &iq_mutex); pthread_mutex_unlock(&iq_mutex); } } } void fl2k_callback(fl2k_data_info_t *data_info) { if (data_info->device_error) { fprintf(stderr, "Device error, exiting.\n"); do_exit = 1; pthread_mutex_lock(&iq_mutex); pthread_cond_signal(&iq_cond); pthread_mutex_unlock(&iq_mutex); } pthread_cond_signal(&cb_cond); data_info->sampletype_signed = 1; data_info->r_buf = (char *)itxbuf; data_info->g_buf = (char *)qtxbuf; } int main(int argc, char **argv) { int r, opt; uint32_t buf_num = 0; int dev_index = 0; pthread_attr_t attr; char *filename = NULL; int option_index = 0; int input_freq_specified = 0; enum inputType_E input_type = INP_REAL; #ifndef _WIN32 struct sigaction sigact, sigign; #endif static struct option long_options[] = { {0, 0, 0, 0} }; while (1) { opt = getopt_long(argc, argv, "ewd:c:i:s:t:", long_options, &option_index); /* end of options reached */ if (opt == -1) break; switch (opt) { case 0: break; case 'd': dev_index = (uint32_t)atoi(optarg); break; case 'c': base_freq = (uint32_t)atof(optarg); break; case 'i': input_freq = (uint32_t)atof(optarg); input_freq_specified = 1; break; case 's': samp_rate = (uint32_t)atof(optarg); break; case 't': /* type */ if(strcasecmp(optarg, "complex") && strcasecmp(optarg, "real")){ fprintf(stderr, "Unknown parameter to -t : %s", optarg); exit(1); } input_type = strcasecmp(optarg, "complex") == 0 ? INP_COMPLEX : INP_REAL; break; case 'w': swap_iq = 1; break; case 'e': ignore_eof = 1; break; default: usage(); break; } } if (argc <= optind) { usage(); } else { filename = argv[optind]; } if (dev_index < 0) { exit(1); } if (strcmp(filename, "-") == 0) { /* Read samples from stdin */ file = stdin; #ifdef _WIN32 _setmode(_fileno(stdin), _O_BINARY); #endif } else { file = fopen(filename, "rb"); if (!file) { fprintf(stderr, "Failed to open %s\n", filename); return -ENOENT; } } /* allocate I buffer */ buf1 = malloc(FL2K_BUF_LEN); buf2 = malloc(FL2K_BUF_LEN); if (!buf1 || !buf2) { fprintf(stderr, "malloc error!\n"); exit(1); } iambuf = buf1; itxbuf = buf2; /* allocate Q buffer */ buf1 = malloc(FL2K_BUF_LEN); buf2 = malloc(FL2K_BUF_LEN); if (!buf1 || !buf2) { fprintf(stderr, "malloc error!\n"); exit(1); } qambuf = buf1; qtxbuf = buf2; /* Baseband buffer */ slopebuf = malloc(BUFFER_SAMPLES * sizeof(float complex)); ampbuf = malloc(BUFFER_SAMPLES * sizeof(float complex)); pdbuf = malloc(BUFFER_SAMPLES * sizeof(long int)); pdslopebuf = malloc(BUFFER_SAMPLES * sizeof(long int)); readpos = 0; writepos = 1; fprintf(stdout, "Samplerate: %3.2f MHz\n", (float)samp_rate/1000000); fprintf(stdout, "Center frequency: %5.0f kHz\n", (float)base_freq/1000); if(swap_iq) fprintf(stdout, "Spectral inversion active.\n"); if(ignore_eof) fprintf(stdout, "Ignoring EOF.\n"); pthread_mutex_init(&cb_mutex, NULL); pthread_mutex_init(&iq_mutex, NULL); pthread_cond_init(&cb_cond, NULL); pthread_cond_init(&iq_cond, NULL); pthread_attr_init(&attr); fl2k_open(&dev, (uint32_t)dev_index); if (NULL == dev) { fprintf(stderr, "Failed to open fl2k device #%d.\n", dev_index); goto out; } r = pthread_create(&iq_thread, &attr, iq_worker, NULL); if (r < 0) { fprintf(stderr, "Error spawning IQ worker thread!\n"); goto out; } pthread_attr_destroy(&attr); r = fl2k_start_tx(dev, fl2k_callback, NULL, 0); /* Set the sample rate */ r = fl2k_set_sample_rate(dev, samp_rate); if (r < 0) fprintf(stderr, "WARNING: Failed to set sample rate. %d\n", r); /* read back actual frequency */ samp_rate = fl2k_get_sample_rate(dev); /* Calculate needed constants */ rf_to_baseband_sample_ratio = samp_rate / input_freq; #ifndef _WIN32 sigact.sa_handler = sighandler; sigemptyset(&sigact.sa_mask); sigact.sa_flags = 0; sigign.sa_handler = SIG_IGN; sigaction(SIGINT, &sigact, NULL); sigaction(SIGTERM, &sigact, NULL); sigaction(SIGQUIT, &sigact, NULL); sigaction(SIGPIPE, &sigign, NULL); #else SetConsoleCtrlHandler( (PHANDLER_ROUTINE) sighandler, TRUE ); #endif ampliphase_modulator(input_type); out: fl2k_close(dev); if (file != stdin) fclose(file); free(ampbuf); free(slopebuf); free(buf1); free(buf2); return 0; }