radiosonde_auto_rx/iq/shift_IQ.c

178 wiersze
4.6 KiB
C
Czysty Zwykły widok Historia

2015-11-16 16:54:41 +00:00
/*
2015-12-02 13:20:38 +00:00
* ./shift IQ.wav [freq] [> IQ_translated.wav] ## translate
* ./shift IQ.wav - [> IQ_swap.wav] ## swap IQ
2015-11-16 16:54:41 +00:00
*/
/*
2015-12-03 19:44:27 +00:00
sdr#<rev1381: sdrsharpIQ.wav = conj(rtlsdrIQ.wav) -> rtlsdrIQ.wav: swapIQ
2015-11-16 16:54:41 +00:00
*/
/*
time domain frequency domain
w(t) = exp(t*i*2pi*f0)*z(t) <-> W(f) = Z(f-f0)
z = x + iy
w = u + iv = exp(ia)*z
u = x*cos(a) - y*sin(a);
v = x*sin(a) + y*cos(a);
*/
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <complex.h>
#define PI (3.1415926535897932384626433832795)
int sample_rate, bits_sample;
int findstr(char *buff, char *str, int pos) {
int i;
for (i = 0; i < 4; i++) {
if (buff[(pos+i)%4] != str[i]) break;
}
return i;
}
int read_IQwavheader(FILE *fp, int *sample_rate, int *bits_sample, FILE *fout) {
int channels = 0;
char txt[4+1] = "\0\0\0\0";
unsigned char dat[4];
int byte, p=0;
if (fread(txt, 1, 4, fp) < 4) return -1; else fwrite(txt, 1, 4, fout);
if (strncmp(txt, "RIFF", 4)) return -1;
// pos_WAVE = 8L
if (fread(txt, 1, 4, fp) < 4) return -1; else fwrite(txt, 1, 4, fout);
if (fread(txt, 1, 4, fp) < 4) return -1; else fwrite(txt, 1, 4, fout);
if (strncmp(txt, "WAVE", 4)) return -1;
// pos_fmt = 12L
for ( ; ; ) {
if ( (byte=fgetc(fp)) == EOF ) return -1; else fprintf(fout, "%c", byte & 0xFF);
txt[p % 4] = byte;
p++; if (p==4) p=0;
if (findstr(txt, "fmt ", p) == 4) break;
}
if (fread(dat, 1, 4, fp) < 4) return -1; else fwrite(dat, 1, 4, fout);
if (fread(dat, 1, 2, fp) < 2) return -1; else fwrite(dat, 1, 2, fout);
if (fread(dat, 1, 2, fp) < 2) return -1; else fwrite(dat, 1, 2, fout);
channels = dat[0] + (dat[1] << 8);
if (fread(dat, 1, 4, fp) < 4) return -1; else fwrite(dat, 1, 4, fout);
*sample_rate = dat[0] | (dat[1] << 8) | (dat[2] << 16) | (dat[3] << 24); //memcpy(&sr, dat, 4);
if (fread(dat, 1, 4, fp) < 4) return -1; else fwrite(dat, 1, 4, fout);
if (fread(dat, 1, 2, fp) < 2) return -1; else fwrite(dat, 1, 2, fout);
if (fread(dat, 1, 2, fp) < 2) return -1; else fwrite(dat, 1, 2, fout);
*bits_sample = dat[0] + (dat[1] << 8);
// pos_dat = 36L + info
for ( ; ; ) {
if ( (byte=fgetc(fp)) == EOF ) return -1; else fprintf(fout, "%c", byte & 0xFF);
txt[p % 4] = byte;
p++; if (p==4) p=0;
if (findstr(txt, "data", p) == 4) break;
}
if (fread(dat, 1, 4, fp) < 4) return -1; else fwrite(dat, 1, 4, fout);
fprintf(stderr, "sample_rate: %d\n", *sample_rate);
fprintf(stderr, "bits : %d\n", *bits_sample);
fprintf(stderr, "channels : %d\n", channels);
if ((*bits_sample != 8) && (*bits_sample != 16)) return -2;
if (channels != 2) return -3;
return 0;
}
int read_csample(FILE *fp, double complex *z) {
2016-01-07 15:33:55 +00:00
short x = 0, y = 0;
2015-11-16 16:54:41 +00:00
if (fread( &x, bits_sample/8, 1, fp) != 1) return EOF;
if (fread( &y, bits_sample/8, 1, fp) != 1) return EOF;
*z = x + I*y;
if (bits_sample == 8) {
*z -= 128 + I*128;
}
return 0;
}
int write_csample(FILE *fp, double complex w) {
int u, v;
if (bits_sample == 8) {
w += 128 + I*128;
}
u = creal(w);
v = cimag(w);
// 16 bit (short) -> (int)
fwrite( &u, bits_sample/8, 1, fp); // + 0000 .. 7FFF -> 0000 0000 .. 0000 7FFF
fwrite( &v, bits_sample/8, 1, fp); // - 8000 .. FFFF -> FFFF 8000 .. FFFF FFFF
return 0;
}
int main(int argc, char *argv[]) {
FILE *fp, *fout;
int sample = 0;
double t, f = 0;
double complex z, w;
2015-12-02 13:20:38 +00:00
int swap = 0;
2015-11-16 16:54:41 +00:00
if (argv[1]) {
fp = fopen(argv[1], "rb");
if (fp == NULL) return -1;
fout = stdout; //fopen("tmp_out.wav", "wb");
2015-12-02 13:20:38 +00:00
if (read_IQwavheader(fp, &sample_rate, &bits_sample, fout) != 0) {
fprintf(stderr, "error: wav header\n");
return -1;
}
2015-11-16 16:54:41 +00:00
2015-12-02 13:20:38 +00:00
if (argv[2]) {
if (strcmp(argv[2], "-") == 0) swap = 1;
else f = (double)atoi(argv[2]); // f0
}
2015-11-16 16:54:41 +00:00
sample = 0;
for ( ; ; ) {
if (read_csample(fp, &z) == EOF) break;
2015-12-02 13:20:38 +00:00
if (swap) { // swap IQ
2015-12-02 13:58:04 +00:00
w = cimag(z) + I*creal(z); // = cexp(I*PI/2) * conj(z) = I * conj(z);
// W(f) = Z(-f); fuer Spiegelung reicht: w = conj(z)
2015-12-02 13:20:38 +00:00
}
else { // translate
t = (double)sample / sample_rate;
w = cexp(t*2*PI*f*I) * z;
}
2015-11-16 16:54:41 +00:00
write_csample(fout, w);
sample++;
}
fclose(fp);
fclose(fout);
}
return 0;
}