diff --git a/scan/plot_fft_pow.py b/scan/plot_fft_pow.py new file mode 100644 index 0000000..b405c8a --- /dev/null +++ b/scan/plot_fft_pow.py @@ -0,0 +1,32 @@ + +import numpy as np +import matplotlib.pyplot as plt + + +# { b: blue , g: green , r: red , c: cyan , m: magenta , y: yellow , k: black , w: white } + +data = np.genfromtxt( 'db.txt', delimiter=';', names=['1','2','3','4','5','6','7','8', '9'] ) + +fq = data['1'] +freq = data['2'] +db = data['3'] +intdb = data['4'] +avg = data['5'] +dev = data['6'] +peak = data['7'] +peak2 = data['8'] +mag = data['9'] + + +ax1 = plt.subplot(211) +ax1.plot( fq, db, color='b', linewidth=0.2, label='fft') +ax1.plot( fq, intdb, color='g', label='fft') + +ax2 = plt.subplot(212) +ax2.plot( fq, avg, color='b', label='fft') +ax2.plot( fq, peak, color='r', label='fft') +ax2.plot( fq, peak2, color='c', label='fft') +ax2.plot( fq, mag, color='y', label='fft') + +plt.show() + diff --git a/scan/scan_fft_pow.c b/scan/scan_fft_pow.c new file mode 100644 index 0000000..24c0e5a --- /dev/null +++ b/scan/scan_fft_pow.c @@ -0,0 +1,617 @@ + +#include +#include +#include +#include +#include + + +typedef unsigned char ui8_t; +typedef unsigned short ui16_t; +typedef unsigned int ui32_t; +typedef short i16_t; +typedef int i32_t; + + +static int option_verbose = 0, // ausfuehrliche Anzeige + option_silent = 0, + wavloaded = 0; +static int wav_channel = 0; // audio channel: left + + +static int sample_rate = 0, bits_sample = 0, channels = 0; +static int wav_ch = 0; // 0: links bzw. mono; 1: rechts + +static unsigned int sample; + +float complex *buffer = NULL; + +void *bufIQ; + +/* ------------------------------------------------------------------------------------ */ + +typedef struct { + int sr; // sample_rate + int LOG2N; + int N; + int N2; + float *xn; + float complex *ew; + float complex *Fm; + float complex *X; + float complex *Z; + float complex *cx; + float complex *win; // float real +} dft_t; + + +dft_t DFT; + +static float *db, *sum_db, *intdb; +static float *avg, *peak; + +static void raw_dft(dft_t *dft, float complex *Z) { + int s, l, l2, i, j, k; + float complex w1, w2, T; + + j = 1; + for (i = 1; i < dft->N; i++) { + if (i < j) { + T = Z[j-1]; + Z[j-1] = Z[i-1]; + Z[i-1] = T; + } + k = dft->N/2; + while (k < j) { + j = j - k; + k = k/2; + } + j = j + k; + } + + for (s = 0; s < dft->LOG2N; s++) { + l2 = 1 << s; + l = l2 << 1; + w1 = (float complex)1.0; + w2 = dft->ew[s]; // cexp(-I*M_PI/(float)l2) + for (j = 1; j <= l2; j++) { + for (i = j; i <= dft->N; i += l) { + k = i + l2; + T = Z[k-1] * w1; + Z[k-1] = Z[i-1] - T; + Z[i-1] = Z[i-1] + T; + } + w1 = w1 * w2; + } + } +} + +static void cdft(dft_t *dft, float complex *z, float complex *Z) { + int i; + for (i = 0; i < dft->N; i++) Z[i] = z[i]; + raw_dft(dft, Z); +} + +static void rdft(dft_t *dft, float *x, float complex *Z) { + int i; + for (i = 0; i < dft->N; i++) Z[i] = (float complex)x[i]; + raw_dft(dft, Z); +} + +static void Nidft(dft_t *dft, float complex *Z, float complex *z) { + int i; + for (i = 0; i < dft->N; i++) z[i] = conj(Z[i]); + raw_dft(dft, z); + // idft(): + // for (i = 0; i < dft->N; i++) z[i] = conj(z[i])/(float)dft->N; // hier: z reell +} + +static float bin2freq(dft_t *dft, int k) { + float fq = k / (float)dft->N; + if ( fq >= 0.5) fq -= 1.0; + return fq*dft->sr; +} + +static float bin2fq(dft_t *dft, int k) { + float fq = k / (float)dft->N; + if ( fq >= 0.5) fq -= 1.0; + return fq; +} + +static float freq2bin(dft_t *dft, int f) { + return f/(float)dft->sr * dft->N; +} + + +/* ------------------------------------------------------------------------------------ */ + +static int dft_window(dft_t *dft, int w) { + int n; + + if (w < 0 || w > 3) return -1; + + for (n = 0; n < dft->N2; n++) { + switch (w) + { + case 0: // (boxcar) + dft->win[n] = 1.0; + break; + case 1: // Hann + dft->win[n] = 0.5 * ( 1.0 - cos(2*M_PI*n/(float)(dft->N2-1)) ); + break ; + case 2: // Hamming + dft->win[n] = 25/46.0 + (1.0 - 25/46.0)*cos(2*M_PI*n / (float)(dft->N2-1)); + break ; + case 3: // Blackmann + dft->win[n] = 7938/18608.0 + - 9240/18608.0*cos(2*M_PI*n / (float)(dft->N2-1)) + + 1430/18608.0*cos(4*M_PI*n / (float)(dft->N2-1)); + break ; + } + } + while (n < dft->N) dft->win[n++] = 0.0; + + return 0; +} + +//static double ilog102 = 0.434294482/2.0; // log(10)/2 +static void db_power(dft_t *dft, float complex Z[], float db[]) { // iq-samples/V [-1..1] + int i; // dBw = 2*dBv, P=c*U*U + for (i = 0; i < dft->N; i++) { // dBw = 2*10*log10(V/V0) + //db[i] = 10.0*log10(creal(Z[i])*creal(Z[i])+cimag(Z[i])*cimag(Z[i])+1e-20); + db[i] = 20.0 * log10(cabs(Z[i])/dft->N2+1e-20); // 20log10(Z/N) + } +} + + +static int init_dft(dft_t *dft) { + int i, k, n; + float normM = 0; + + bufIQ = calloc(2*(dft->N+2), 2); if (bufIQ == NULL) return -1; + buffer = calloc(dft->N+1, sizeof(float complex)); if (buffer == NULL) return -1; + + dft->xn = calloc(dft->N+1, sizeof(float complex)); if (dft->xn == NULL) return -1; + dft->Z = calloc(dft->N+1, sizeof(float complex)); if (dft->Z == NULL) return -1; + dft->ew = calloc(dft->LOG2N, sizeof(float complex)); if (dft->ew == NULL) return -1; + + dft->win = calloc(dft->N+1, sizeof(float complex)); if (dft->win == NULL) return -1; + dft->N2 = dft->N; + dft_window(dft, 1); + + normM = 0; + for (i = 0; i < dft->N2; i++) normM += dft->win[i]*dft->win[i]; + //normM = sqrt(normM); + + for (n = 0; n < dft->LOG2N; n++) { + k = 1 << n; + dft->ew[n] = cexp(-I*M_PI/(double)k); + } + + db = calloc(dft->N+1, sizeof(float)); if (db == NULL) return -1; + sum_db = calloc(dft->N+1, sizeof(float)); if (sum_db == NULL) return -1; + intdb = calloc(dft->N+1, sizeof(float)); if (intdb == NULL) return -1; + avg = calloc(dft->N+1, sizeof(float)); if (avg == NULL) return -1; + peak = calloc(dft->N+1, sizeof(float)); if (peak == NULL) return -1; + + return 0; +} + + +static void end_dft(dft_t *dft) { + if (bufIQ) { free(bufIQ); bufIQ = NULL; } + if (buffer) { free(buffer); buffer = NULL; } + if (dft->xn) { free(dft->xn); dft->xn = NULL; } + if (dft->Z) { free(dft->Z); dft->Z = NULL; } + if (dft->ew) { free(dft->ew); dft->ew = NULL; } + if (dft->win) { free(dft->win); dft->win = NULL; } + if (db) { free(db); db = NULL; } + if (sum_db) { free(sum_db); sum_db = NULL; } + if (intdb) { free(intdb); intdb = NULL; } + if (avg) { free(avg); avg = NULL; } + if (peak) { free(peak); peak = NULL; } +} + + +/* ------------------------------------------------------------------------------------ */ + +static 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; +} + +static int read_wav_header(FILE *fp, int wav_channel) { + 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; + if (strncmp(txt, "RIFF", 4)) return -1; + if (fread(txt, 1, 4, fp) < 4) return -1; + // pos_WAVE = 8L + if (fread(txt, 1, 4, fp) < 4) return -1; + if (strncmp(txt, "WAVE", 4)) return -1; + // pos_fmt = 12L + for ( ; ; ) { + if ( (byte=fgetc(fp)) == EOF ) return -1; + 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; + if (fread(dat, 1, 2, fp) < 2) return -1; + + if (fread(dat, 1, 2, fp) < 2) return -1; + channels = dat[0] + (dat[1] << 8); + + if (fread(dat, 1, 4, fp) < 4) return -1; + memcpy(&sample_rate, dat, 4); //sample_rate = dat[0]|(dat[1]<<8)|(dat[2]<<16)|(dat[3]<<24); + + if (fread(dat, 1, 4, fp) < 4) return -1; + if (fread(dat, 1, 2, fp) < 2) return -1; + //byte = dat[0] + (dat[1] << 8); + + if (fread(dat, 1, 2, fp) < 2) return -1; + bits_sample = dat[0] + (dat[1] << 8); + + // pos_dat = 36L + info + for ( ; ; ) { + if ( (byte=fgetc(fp)) == EOF ) return -1; + 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; + + + fprintf(stderr, "sample_rate: %d\n", sample_rate); + fprintf(stderr, "bits : %d\n", bits_sample); + fprintf(stderr, "channels : %d\n", channels); + + if (wav_channel >= 0 && wav_channel < channels) wav_ch = wav_channel; + else wav_ch = 0; + fprintf(stderr, "channel-In : %d\n", wav_ch+1); + + if ((bits_sample != 8) && (bits_sample != 16)) return -1; + + return 0; +} + + +static int read_bufIQ(dft_t *dft, FILE *fp) { + int len; + + len = fread( bufIQ, bits_sample/8, 2*dft->N2, fp); + + if ( len != 2*dft->N2) { + while (len < 2*dft->N2) { + //bufIQ[len] = 0; + len++; + } + return EOF; + } + + return 0; +} + +static int bufIQ2complex(dft_t *dft) { + int i; + float complex z; + unsigned char *buf8; + short *buf16; + + if (bits_sample == 8) { + buf8 = bufIQ; + for (i = 0; i < dft->N2; i++) { + z = buf8[2*i]-128.0 + I*(buf8[2*i+1]-128.0); + z /= 128.0; + buffer[i] = z; + } + } + else { // bits_sample == 16 + buf16 = bufIQ; + for (i = 0; i < dft->N2; i++) { + z = buf16[2*i] + I*buf16[2*i+1]; + z /= 128.0*256.0; + buffer[i] = z; + } + } + + return 0; +} + +static int f32read_sample(FILE *fp, float *s) { + int i; + short b = 0; + + for (i = 0; i < channels; i++) { + + if (fread( &b, bits_sample/8, 1, fp) != 1) return EOF; + + if (i == wav_ch) { // i = 0: links bzw. mono + //if (bits_sample == 8) sint = b-128; // 8bit: 00..FF, centerpoint 0x80=128 + //if (bits_sample == 16) sint = (short)b; + + if (bits_sample == 8) { b -= 128; } + *s = b/128.0; + if (bits_sample == 16) { *s /= 256.0; } + } + } + + return 0; +} + +/* ------------------------------------------------------------------------------------ */ + + +int main(int argc, char **argv) { + + FILE *OUT = stderr; + FILE *fpout = stdout; + FILE *fp = NULL; + char *prgnam = NULL; + char *filename = NULL; + + int mn = 0; // 0: N = M + + int j, n, k; + float tl = 4.0; + + float dx; + int dn; + //float avg = 0.0; + float dev = 0.0; + float sympeak = 0.0; + float sympeak2 = 0.0; + float globmin = 0.0; + float globavg = 0.0; + + +#ifdef CYGWIN + _setmode(fileno(stdin), _O_BINARY); // _setmode(_fileno(stdin), _O_BINARY); +#endif + setbuf(stdout, NULL); + + prgnam = argv[0]; + ++argv; + while ((*argv) && (!wavloaded)) { + if ( (strcmp(*argv, "-h") == 0) || (strcmp(*argv, "--help") == 0) ) { + fprintf(stderr, "%s [options] audio.wav\n", prgnam); + fprintf(stderr, " options:\n"); + //fprintf(stderr, " -v, --verbose\n"); + return 0; + } + else if ( (strcmp(*argv, "-v") == 0) || (strcmp(*argv, "--verbose") == 0) ) { + option_verbose = 1; + } + else if ( (strcmp(*argv, "-s") == 0) || (strcmp(*argv, "--silent") == 0) ) { + option_silent = 1; + } + else if ( (strcmp(*argv, "-t") == 0) || (strcmp(*argv, "--time") == 0) ) { + ++argv; + if (*argv) tl = atof(*argv); + else return -1; + } + else { + if (strcmp(*argv, "-") == 0) { + if (argv[1] == NULL) return -1; else sample_rate = atoi(argv[1]); + if (argv[2] == NULL) return -1; else bits_sample = atoi(argv[2]); + channels = 2; + if (argv[3] == NULL) fp = stdin; + else { + fp = fopen(argv[3], "rb"); + if (fp == NULL) { + fprintf(stderr, "%s konnte nicht geoeffnet werden\n", *argv); + return -1; + } + } + wavloaded = 2; + } + else { + fp = fopen(*argv, "rb"); + if (fp == NULL) { + fprintf(stderr, "%s konnte nicht geoeffnet werden\n", *argv); + return -1; + } + wavloaded = 1; + } + } + ++argv; + } + if (!wavloaded) fp = stdin; + + + if (wavloaded < 2) { + j = read_wav_header(fp, wav_channel); + if ( j < 0 ) { + fclose(fp); + fprintf(stderr, "error: wav header\n"); + return -1; + } + } + + DFT.sr = sample_rate; + DFT.LOG2N = 14; // 2^12=4096: 300-400Hz bins, 2^14=16384: 75-100 Hz + mn = 0; + DFT.N2 = 1 << DFT.LOG2N; + if (DFT.N2 > DFT.sr/2) { + DFT.LOG2N = 0; + while ( (1 << (DFT.LOG2N+1)) < DFT.sr/2 ) DFT.LOG2N++; + DFT.N2 = 1 << DFT.LOG2N; + } + DFT.N = DFT.N2 << mn; + DFT.LOG2N += mn; + + init_dft(&DFT); + + if (option_verbose) fprintf(stderr, "M: %d\n", DFT.N2); + + + //memset(sum_db, 0, N*sizeof(float)); // calloc() + + sample = 0; + n = 0; + while ( read_bufIQ(&DFT, fp) != EOF ) { + + bufIQ2complex(&DFT); + sample += DFT.N2; + + if (sample % DFT.N2 == 0) + { + double complex dc = 0; + for (j = 0; j < DFT.N; j++) { + DFT.Z[j] = buffer[j % DFT.N]; + dc += DFT.Z[j]; + } + dc /= 0.99*DFT.N; + //dc = 0; + + for (j = 0; j < DFT.N; j++) { + DFT.Z[j] = buffer[j % DFT.N] - dc; + } +///* + for (j = 0; j < DFT.N2; j++) { + DFT.Z[j] *= DFT.win[j]; + } + while (j < DFT.N) DFT.Z[j++] = 0.0; +//*/ + raw_dft(&DFT, DFT.Z); + db_power(&DFT, DFT.Z, db); + + for (j = 0; j < DFT.N; j++) sum_db[j] += db[j]; + + n++; + + if (sample > tl*DFT.sr) break; + } + } + if (option_verbose) fprintf(stderr, "n=%d\n", n); + + if (option_verbose == 0) { + OUT = stdout; + fpout = fopen("db.txt", "wb"); + if (fpout == NULL) return -1; + } else { + OUT = stderr; + fpout = stdout; + } + + + globmin = 0.0; + globavg = 0.0; + float db_spike3 = 10.0; + int spike_wl3 = 3; //freq2bin(&DFT, 200); // 3 // 200 Hz + int spike_wl5 = 5; //freq2bin(&DFT, 200); // 3 // 200 Hz + float db_spike1 = 15.0; + int spike_wl1 = 1; //freq2bin(&DFT, 200); // 3 // 200 Hz + //fprintf(stderr, "sp_wl: %d\n", spike_wl); + + dx = bin2freq(&DFT, 1); + dn = 2*(int)(2400.0/dx)+1; // (odd/symmetric) integration width: 4800+dx Hz + if (option_verbose) fprintf(stderr, "dn = %d\n", dn); + for (j = 0; j < DFT.N; j++) sum_db[j] /= (float)n; + + // dc-spike (N-1,)N,0,1(,2): subtract mean/avg + // spikes in general: + for (j = 0; j < DFT.N; j++) { + if ( sum_db[j] - sum_db[(j-spike_wl5+DFT.N)%DFT.N] > db_spike3 + && sum_db[j] - sum_db[(j-spike_wl3+DFT.N)%DFT.N] > db_spike3 + && sum_db[j] - sum_db[(j+spike_wl3+DFT.N)%DFT.N] > db_spike3 + && sum_db[j] - sum_db[(j+spike_wl5+DFT.N)%DFT.N] > db_spike3 + ) { + sum_db[j] = (sum_db[(j-spike_wl3+DFT.N)%DFT.N]+sum_db[(j+spike_wl3+DFT.N)%DFT.N])/2.0; + } + } + + for (j = 0; j < DFT.N; j++) { + float sum = 0.0; + for (n = j-(dn-1)/2; n <= j+(dn-1)/2; n++) sum += sum_db[(n + DFT.N) % DFT.N]; + sum /= (float)dn; + intdb[j] = sum; + globavg += sum; // <=> sum_db[j]; + if (sum < globmin) globmin = sum; + } + globavg /= (float)DFT.N; + + if (option_verbose) fprintf(stderr, "avg=%.2f\n", globavg); + if (option_verbose) fprintf(stderr, "min=%.2f\n", globmin); + + int dn2 = 2*dn+1; + int dn3 = (int)(4000.0/dx); // (odd/symmetric) integration width: +/-4000 Hz + + int delay = (int)(24000.0/dx); // 16000 + float *bufavg = calloc(delay+1, sizeof(float)); if (bufavg == NULL) return -1; + float *bufpeak = calloc(delay+1, sizeof(float)); if (bufpeak == NULL) return -1; + k = 0; + int mag = 0; + int mag0 = 0; + float max_db_loc = 0.0; + int max_db_idx = 0; + for (j = DFT.N/2; j < DFT.N/2 + DFT.N; j++) { + float x = intdb[(j+delay) % DFT.N]; + float a = 0.0; + + if (1 || option_verbose) fprintf(fpout, "%9.6f ; %9.1f ; %10.4f", bin2fq(&DFT, j % DFT.N), bin2freq(&DFT, j % DFT.N), sum_db[j % DFT.N]); + if (1 || option_verbose) fprintf(fpout, " ; %10.4f", intdb[j % DFT.N]); + + a = 0.0; + for (n = j-(dn2-1)/2; n <= j+(dn2-1)/2; n++) a += intdb[n % DFT.N]; + a /= (float)dn2; + bufavg[k % delay] = a; + + dev = 0.0; + for (n = j-(dn2-1)/2; n <= j+(dn2-1)/2; n++) dev += (intdb[n % DFT.N]-a)*(intdb[n % DFT.N]-a); + dev = sqrt(dev/(float)dn2); + + + sympeak = 0.0; + for (n = 1; n <= dn3; n++) { + sympeak += (sum_db[(j+n) % DFT.N]-globmin)*(sum_db[(j-n + DFT.N) % DFT.N]-globmin); + } + sympeak = sqrt(abs(sympeak)/(float)dn3); // globmin > min + + sympeak2 = 0.0; + for (n = 0; n <= (dn2-1)/2; n++) sympeak2 += (intdb[(j+n) % DFT.N]-globmin)*(intdb[(j-n + DFT.N) % DFT.N]-globmin); + sympeak2 = sqrt(sympeak2/(2.0*dn2)); + + peak[k % delay] = sympeak; + + + if (1 || option_verbose) fprintf(fpout, " ; %10.4f ; %10.4f", a-globmin, dev); + if (1 || option_verbose) fprintf(fpout, " ; %10.4f ; %10.4f", sympeak, sympeak2); + + mag = (sympeak - peak[(k+1)%delay])/ 3.0; // threshold 3.0 + if ( mag < 0 ) mag = 0; + if (1 || option_verbose) fprintf(fpout, " ; %d", mag); + + if (mag0 > 0 && mag == 0) fprintf(OUT, "peak: %+9.6f = %+9.1fHz\n", bin2fq(&DFT, max_db_idx), bin2freq(&DFT, max_db_idx)); + if (mag0 == 0 && mag > 0) { + max_db_loc = sympeak; + max_db_idx = j % DFT.N; + } + if (mag > 0 && sympeak > max_db_loc) { + max_db_loc = sympeak; + max_db_idx = j % DFT.N; + } + + mag0 = mag; + if (1 || option_verbose) fprintf(fpout, "\n"); + + k++; + } + + if (option_verbose) fprintf(stderr, "bin = %.2f Hz\n", bin2freq(&DFT, 1)); + + if (bufavg) free(bufavg); + if (bufpeak) free(bufpeak); + + end_dft(&DFT); + fclose(fp); + + + return 0; +} + diff --git a/scan/scan_multi.sh b/scan/scan_multi.sh new file mode 100755 index 0000000..29200a2 --- /dev/null +++ b/scan/scan_multi.sh @@ -0,0 +1,16 @@ + +#mkfifo /tmp/sdr.0 +#mkfifo /tmp/sdr.1 +#mkfifo /tmp/sdr.2 + +#mkfifo log.0 +#mkfifo log.1 +#mkfifo log.2 + +xfce4-terminal --geometry=64x61+470+140 --hide-menubar -T term0 -H -e './scan3term.pl' + +xfce4-terminal --geometry=160x20+1000+836 --hide-menubar -T term3 -H -e 'cat log.2' +xfce4-terminal --geometry=160x20+1000+488 --hide-menubar -T term2 -H -e 'cat log.1' +xfce4-terminal --geometry=160x20+1000+140 --hide-menubar -T term1 -H -e 'cat log.0' + + diff --git a/scan/scan_multi_rs.pl b/scan/scan_multi_rs.pl new file mode 100755 index 0000000..8a7f422 --- /dev/null +++ b/scan/scan_multi_rs.pl @@ -0,0 +1,200 @@ +#!/usr/bin/env perl + +use strict; +use warnings; +use 5.010; + + +my $verbose = 1; +my $stderr_null = " 2>/dev/null "; + +my $j; +my $max_terms = 3; +# mkfifo /tmp/sdr.X +# mkfifo log.X , X=0,1,2 +my $log = "log"; +my $sdr = "/tmp/sdr"; + +my $ppm = 45; +my $gain = 40.2; + +my $if_sr = 48e3; +my $decimate = 40; # 20..50 # 32*48000 = 1536000 +my $band_sr = $if_sr * $decimate; # 40*48000 = 1920000 + +my $center_freq = 404620e3; # off-center, rtl_sdr-mirrors ... + +if ($ARGV[0]) { + $band_sr = $ARGV[0]; + if ($ARGV[1]) { + $decimate = $ARGV[1]; + if ($decimate < 20) { $decimate = 20; } + if ($decimate > 50) { $decimate = 50; } + $if_sr = $band_sr / $decimate; + if ($ARGV[2]) { + $center_freq = $ARGV[2]; + } + } +} + + +print "\n"; + +print "[scan ".($center_freq-$band_sr/2).":".($center_freq).":".($center_freq+$band_sr/2)."]\n"; +print "\n"; + +my $str_rtlsdr = "rtl_sdr -g $gain -p ".$ppm." -f ".$center_freq." -s ".$band_sr; +print $str_rtlsdr."\n"; + +my $bps = 8; # rtl-sdr: 8bit +my $b2f_conv; +if ($bps == 16) { $b2f_conv = "convert_s16_f"; } +else { $b2f_conv = "convert_u8_f"; } + + +print "\n"; + + +my $iqraw = "rtlsdr.raw"; +system("$str_rtlsdr -n ".(8*$band_sr)." $iqraw"); + +my $powfile = "peaks.txt"; + +system("./scan_pow - $band_sr $bps $iqraw > $powfile"); +system("python plot_fft_pow.py &"); + + +my @rs_matrix; + +$rs_matrix[2][0] = "dfm"; +$rs_matrix[2][1] = "9600"; # lp-filter bw +$rs_matrix[2][2] = "./dfm09 --ecc --ptu -v"; # --auto ## DFM9 = -DFM + +$rs_matrix[3][0] = "rs41"; +$rs_matrix[3][1] = "9600"; # lp-filter bw +$rs_matrix[3][2] = "./rs41 --ecc2 --crc --ptu -vx"; + +$rs_matrix[4][0] = "rs92"; +$rs_matrix[4][1] = "9600"; +$rs_matrix[4][2] = "./rs92 --ecc --crc -v"; # -e brdc / -a almanac + +$rs_matrix[8][0] = "lms6"; +$rs_matrix[8][1] = "9600"; +$rs_matrix[8][2] = "./lms6 --vit --ecc -v"; + +$rs_matrix[5][0] = "m10"; # scan: carrier offset +$rs_matrix[5][1] = "9600"; +$rs_matrix[5][2] = "./m10 -c -vv"; + +$rs_matrix[6][0] = "imet1ab"; +$rs_matrix[6][1] = "32000"; # > detect-bw +$rs_matrix[6][2] = "./imet1ab_dft -v"; + +$rs_matrix[7][0] = "imet4"; +$rs_matrix[7][1] = "12000"; +$rs_matrix[7][2] = "./imet1rs_dft -r"; + +$rs_matrix[10][0] = ""; + +my $detect = "./dft_detect"; + + +my $fh; +open ($fh, '<', "$powfile") or die "error open '$powfile': $!\n"; + +my @rs_array = (); +my @peakarray = (); +my $num_peaks; + +my $num_lines = 0; +my $line; + +my $ret; +my $filter = 12e3; # 12kHz detect-bw + +while ($line = <$fh>) { + $num_lines += 1; + + if ( ($line =~ /peak: *([+-]?\d*\.\d*) = .*/) ) + { + my $fq = $1; + if ($verbose) { print "[ $fq ] "; } + print $line; # no chomp + my $lp = $filter/($if_sr*2.0); + my $cmd = "cat $iqraw | csdr $b2f_conv | csdr shift_addition_cc ".(-$fq)." $stderr_null| ". + "csdr fir_decimate_cc $decimate 0.005 $stderr_null| csdr bandpass_fir_fft_cc -$lp $lp 0.02 $stderr_null| ". + "csdr fmdemod_quadri_cf | csdr limit_ff | csdr convert_f_s16 | ". + "sox -t raw -r $if_sr -b 16 -e signed - -t wav - 2>/dev/null | $detect -t 10 2>/dev/null"; + #if ($verbose) { print $cmd."\n"; } + system("$cmd"); + $ret = $? >> 8; + if ($ret & 0x80) { + $ret = - (0x100 - $ret); # ($ret & 0x80) = core dump + } + if ($ret) { + push @peakarray, $fq; + push @rs_array, $ret; + print "\n"; + } + } +} + + +$num_peaks = scalar(@peakarray); +print "rs-peaks: ".$num_peaks."\n"; + +my $rs; +my $decoder; + + +# mkfifo /tmp/sdr.X ... +# mkfifo log.X + +my $tee = " > $sdr.0"; + +for ($j = 0; $j < $num_peaks; $j++) { + if ($j < $max_terms) { + if ($j > 0) { + $tee = " | tee $sdr.$j ".$tee; + } + + my $idx = abs($rs_array[$j]); + if ($idx > 1 && $idx < 9) { + $rs = $rs_matrix[$idx][0]; + $filter = $rs_matrix[$idx][1]; + my $inv = ($rs_array[$j] < 0 ? "-i" : ""); + if ($idx == 5 || $idx == 7) { $inv = ""; } # || $idx==8 + $decoder = sprintf("%s %s", $rs_matrix[$idx][2], $inv); + + my $lp = $filter/($if_sr*2.0); + + my $pid = fork(); + die if not defined $pid; + if (not $pid) { + my $rs_str = sprintf("%s_%.0fkHz", $rs, ($center_freq+$band_sr*$peakarray[$j])*1e-3); + if ($verbose) { print "\nrs: <".$rs_str.">\n"; } + my $cmd = "cat $sdr.$j | csdr $b2f_conv | csdr shift_addition_cc ".(-$peakarray[$j])." $stderr_null| ". + "csdr fir_decimate_cc $decimate 0.005 $stderr_null| csdr bandpass_fir_fft_cc -$lp $lp 0.02 $stderr_null| ". + "csdr fmdemod_quadri_cf | csdr limit_ff | csdr convert_f_s16 | ". + "sox -t raw -r $if_sr -b 16 -e signed - -t wav - 2>/dev/null | $decoder 2>/dev/null > $log.$j"; + if ($verbose) { print $cmd."\n"; } + system("$cmd"); + exit; + } + } + } +} + + +sleep(1); +print "\n"; + +if ($num_peaks > 0) +{ + my $cmd = "$str_rtlsdr -n ".(60*$band_sr)." - $tee"; + if ($verbose) { print $cmd."\n"; } + system("$cmd"); # timeout +} + +#wait(); +