diff --git a/Makefile b/Makefile index a90dfe9..3e99959 100644 --- a/Makefile +++ b/Makefile @@ -28,7 +28,7 @@ clean: encode: encode.o mmap_file.o pcm.o wav.o alsa.o yuv.o img.o ppm.o sdl.o -decode: decode.o mmap_file.o pcm.o wav.o alsa.o window.o ddc.o delay.o yuv.o img.o ppm.o sdl.o +decode: decode.o mmap_file.o pcm.o wav.o alsa.o window.o ddc.o buffer.o yuv.o img.o ppm.o sdl.o -debug: debug.o mmap_file.o pcm.o wav.o alsa.o window.o ddc.o delay.o yuv.o img.o ppm.o sdl.o +debug: debug.o mmap_file.o pcm.o wav.o alsa.o window.o ddc.o buffer.o yuv.o img.o ppm.o sdl.o diff --git a/delay.c b/buffer.c similarity index 54% rename from delay.c rename to buffer.c index e5e7cc6..e860467 100644 --- a/delay.c +++ b/buffer.c @@ -5,31 +5,34 @@ To the extent possible under law, the author(s) have dedicated all copyright and You should have received a copy of the CC0 Public Domain Dedication along with this software. If not, see . */ - -#include "delay.h" +#include "buffer.h" #include -float do_delay(struct delay *d, float input) +float *do_buffer(struct buffer *d, float input) { - d->s[d->last] = input; - d->last = (d->last + 1) < d->len ? d->last + 1 : 0; - return d->s[d->last]; + d->s[d->last0] = input; + d->s[d->last1] = input; + d->last0 = (d->last0 - 1) < 0 ? d->len : d->last0 - 1; + d->last1 = (d->last1 - 1) < 0 ? d->len : d->last1 - 1; + int last = d->last0 < d->last1 ? d->last0 : d->last1; + return d->s + last; } -struct delay *alloc_delay(int samples) +struct buffer *alloc_buffer(int samples) { - int len = samples + 1; - struct delay *d = malloc(sizeof(struct delay)); + int len = 2 * samples; + struct buffer *d = malloc(sizeof(struct buffer)); d->s = malloc(sizeof(float) * len); - d->last = 0; + d->last0 = 0; + d->last1 = samples; d->len = len; for (int i = 0; i < len; i++) d->s[i] = 0.0; return d; } -void free_delay(struct delay *delay) +void free_buffer(struct buffer *buffer) { - free(delay->s); - free(delay); + free(buffer->s); + free(buffer); } diff --git a/delay.h b/buffer.h similarity index 72% rename from delay.h rename to buffer.h index 64ea1dd..77d5477 100644 --- a/delay.h +++ b/buffer.h @@ -5,17 +5,17 @@ To the extent possible under law, the author(s) have dedicated all copyright and You should have received a copy of the CC0 Public Domain Dedication along with this software. If not, see . */ - -#ifndef DELAY_H -#define DELAY_H -struct delay { +#ifndef BUFFER_H +#define BUFFER_H +struct buffer { float *s; - int last; + int last0; + int last1; int len; }; -float do_delay(struct delay *, float); -struct delay *alloc_delay(int); -void free_delay(struct delay *); +float *do_buffer(struct buffer *d, float input); +struct buffer *alloc_buffer(int samples); +void free_buffer(struct buffer *buffer); #endif diff --git a/ddc.c b/ddc.c index b895065..848da21 100644 --- a/ddc.c +++ b/ddc.c @@ -5,78 +5,70 @@ To the extent possible under law, the author(s) have dedicated all copyright and You should have received a copy of the CC0 Public Domain Dedication along with this software. If not, see . */ - +#include #include #include -#include #include "window.h" #include "ddc.h" -void do_ddc(struct ddc *ddc, float *input, float complex *output) +void do_ddc(struct ddc *ddc, float *input, complex float *output) { - // this works only for L <= M - for (int k = 0, last = ddc->last, in = 0; k < ddc->L; k++) { - while (ddc->skip < ddc->M) { - ddc->s[ddc->last] = input[in++]; - last = ddc->last; - ddc->last = (ddc->last + 1) < ddc->samples ? ddc->last + 1 : 0; - ddc->skip += ddc->L; - } + int N = ddc->N; + int M = ddc->M; + int L = ddc->L; + int offset = 0; + for (int k = 0; k < L; k++) { + float *x = input + (((L * M - 1) - offset) / L); + complex float *b = ddc->b + (N * (offset % L)); + offset += M; - ddc->skip %= ddc->M; - - float complex sum = 0.0; - for (int i = ddc->offset; i < ddc->taps; i += ddc->L) { - sum += ddc->b[i] * ddc->s[last]; - last += last ? - 1 : ddc->samples - 1; - } - - ddc->offset = (ddc->offset + ddc->M) % ddc->L; + complex float sum = 0.0; + for (int i = 0; i < N; i++) + sum += b[i] * x[i]; output[k] = ddc->osc * sum; ddc->osc *= ddc->d; -// ddc->osc /= cabsf(ddc->osc); // not really needed + // ddc->osc /= cabsf(ddc->osc); // not really needed } } -struct ddc *alloc_ddc(float freq, float bw, float step, int taps, int L, int M, float (*window)(float, float, float), float a) +struct ddc *alloc_ddc(int L, int M, float carrier, float bw, float rate, int taps, float (*window)(float, float, float), float a) { - float lstep = step / (float)L; - float ostep = step * (float)M / (float)L; + float lstep = 1.0 / ((float)L * rate); + float ostep = (float)M * lstep; struct ddc *ddc = malloc(sizeof(struct ddc)); - ddc->taps = taps; - ddc->samples = (taps + L - 1) / L; - ddc->b = malloc(sizeof(float complex) * ddc->taps); - ddc->s = malloc(sizeof(float) * ddc->samples); + ddc->N = (taps + L - 1) / L; + ddc->b = malloc(sizeof(complex float) * ddc->N * L); ddc->osc = I; - ddc->d = cexpf(-I * 2.0 * M_PI * freq * ostep); - ddc->offset = (M - 1) % L; - ddc->last = 0; - ddc->skip = 0; + ddc->d = cexpf(-I * 2.0 * M_PI * carrier * ostep); ddc->L = L; ddc->M = M; - for (int i = 0; i < ddc->samples; i++) - ddc->s[i] = 0.0; + complex float *b = malloc(sizeof(complex float) * taps); float sum = 0.0; - for (int i = 0; i < ddc->taps; i++) { - float N = (float)ddc->taps; - float n = (float)i; + for (int i = 0; i < taps; i++) { + float N = taps; + float n = i; float x = n - (N - 1.0) / 2.0; - float l = 2.0 * M_PI * bw * lstep; - float w = window(n, ddc->taps, a); - float h = 0.0 == x ? l / M_PI : sinf(l * x) / (x * M_PI); - float b = w * h; - sum += b; - complex float o = cexpf(I * 2.0 * M_PI * freq * lstep * n); - ddc->b[i] = b * o * (float)L; + float l = 2.0 * bw * lstep; + float h = l * sinc(l * x); + float w = window(n, N, a); + float t = w * h; + sum += t; + complex float o = cexpf(I * 2.0 * M_PI * carrier * lstep * n); + b[i] = t * o; } - for (int i = 0; i < ddc->taps; i++) - ddc->b[i] /= sum; + for (int i = 0; i < taps; i++) + b[i] /= sum; + for (int i = 0; i < ddc->N * L; i++) + ddc->b[i] = 0.0; + for (int i = 0; i < L; i++) + for (int j = i, k = 0; j < taps; j += L, k++) + ddc->b[i * ddc->N + k] = b[j]; + free(b); return ddc; } void free_ddc(struct ddc *ddc) { free(ddc->b); - free(ddc->s); free(ddc); } diff --git a/ddc.h b/ddc.h index e539f60..c3e187d 100644 --- a/ddc.h +++ b/ddc.h @@ -5,28 +5,22 @@ To the extent possible under law, the author(s) have dedicated all copyright and You should have received a copy of the CC0 Public Domain Dedication along with this software. If not, see . */ - #ifndef DDC_H #define DDC_H #include "window.h" struct ddc { - float complex *b; - float *s; - float complex osc; - float complex d; - int offset; - int skip; - int last; - int taps; - int samples; + complex float *b; + complex float osc; + complex float d; + int N; int L; int M; }; -void do_ddc(struct ddc *, float *, float complex *); -struct ddc *alloc_ddc(float, float, float, int, int, int, float (*)(float, float, float), float); -void free_ddc(struct ddc *); +void do_ddc(struct ddc *ddc, float *input, complex float *output); +struct ddc *alloc_ddc(int L, int M, float carrier, float bw, float rate, int taps, float (*window)(float, float, float), float a); +void free_ddc(struct ddc *ddc); #endif diff --git a/debug.c b/debug.c index 37beb5d..4c3fb2b 100644 --- a/debug.c +++ b/debug.c @@ -15,7 +15,7 @@ You should have received a copy of the CC0 Public Domain Dedication along with t #include "mmap_file.h" #include "pcm.h" #include "ddc.h" -#include "delay.h" +#include "buffer.h" #include "yuv.h" #include "utils.h" #include "img.h" @@ -47,7 +47,6 @@ int main(int argc, char **argv) if (channels > 1) fprintf(stderr, "using first of %d channels\n", channels); - const float step = 1.0 / rate; float complex cnt_last = -I; float complex dat_last = -I; @@ -105,18 +104,19 @@ int main(int argc, char **argv) float drate = rate * (float)factor_L / (float)factor_M; float dstep = 1.0 / drate; fprintf(stderr, "using factor of %ld/%ld, working at %.2fhz\n", factor_L, factor_M, drate); - float *cnt_amp = malloc(sizeof(float) * factor_M); - float *dat_amp = malloc(sizeof(float) * factor_M); float complex *cnt_q = malloc(sizeof(float complex) * factor_L); float complex *dat_q = malloc(sizeof(float complex) * factor_L); // same factor to keep life simple and have accurate horizontal sync - struct ddc *cnt_ddc = alloc_ddc(1200.0, 200.0, step, cnt_taps, factor_L, factor_M, kaiser, 2.0); - struct ddc *dat_ddc = alloc_ddc(1900.0, 800.0, step, dat_taps, factor_L, factor_M, kaiser, 2.0); + struct ddc *cnt_ddc = alloc_ddc(factor_L, factor_M, 1200.0, 200.0, rate, cnt_taps, kaiser, 2.0); + struct ddc *dat_ddc = alloc_ddc(factor_L, factor_M, 1900.0, 800.0, rate, dat_taps, kaiser, 2.0); // delay input by phase shift of other filter to synchronize outputs - struct delay *cnt_delay = alloc_delay((dat_taps - 1) / (2 * factor_L)); - struct delay *dat_delay = alloc_delay((cnt_taps - 1) / (2 * factor_L)); + int cnt_delay = (dat_taps - 1) / (2 * factor_L); + int dat_delay = (cnt_taps - 1) / (2 * factor_L); - short *buff = (short *)malloc(sizeof(short) * channels * factor_M); + short *pcm_buff = (short *)malloc(sizeof(short) * channels * factor_M); + + // 0.1 second history + enough room for delay and taps + struct buffer *buffer = alloc_buffer(0.1 * rate + 2 * fmaxf(cnt_delay, dat_delay)); const float sync_porch_len = 0.003; const float porch_len = 0.0015; (void)porch_len; @@ -143,15 +143,14 @@ int main(int argc, char **argv) for (int out = factor_L;; out++, hor_ticks++, cal_ticks++, vis_ticks++) { if (out >= factor_L) { out = 0; - if (!read_pcm(pcm, buff, factor_M)) + if (!read_pcm(pcm, pcm_buff, factor_M)) break; - for (int j = 0; j < factor_M; j++) { - float amp = (float)buff[j * channels] / 32767.0; - cnt_amp[j] = do_delay(cnt_delay, amp); - dat_amp[j] = do_delay(dat_delay, amp); - } - do_ddc(cnt_ddc, cnt_amp, cnt_q); - do_ddc(dat_ddc, dat_amp, dat_q); + float *buff = 0; + for (int j = 0; j < factor_M; j++) + buff = do_buffer(buffer, (float)pcm_buff[j * channels] / 32767.0); + + do_ddc(cnt_ddc, buff + cnt_delay, cnt_q); + do_ddc(dat_ddc, buff + dat_delay, dat_q); } float cnt_freq = fclampf(1200.0 + cargf(cnt_q[out] * conjf(cnt_last)) / (2.0 * M_PI * dstep), 1100.0, 1300.0); @@ -417,8 +416,8 @@ int main(int argc, char **argv) free_ddc(cnt_ddc); free_ddc(dat_ddc); - free(cnt_amp); - free(dat_amp); + free_buffer(buffer); + free(pcm_buff); return 0; } diff --git a/decode.c b/decode.c index 62b8798..33c6848 100644 --- a/decode.c +++ b/decode.c @@ -14,7 +14,7 @@ You should have received a copy of the CC0 Public Domain Dedication along with t #include #include "pcm.h" #include "ddc.h" -#include "delay.h" +#include "buffer.h" #include "yuv.h" #include "utils.h" #include "img.h" @@ -289,23 +289,20 @@ int demodulate(struct pcm *pcm, float *cnt_freq, float *dat_freq, float *drate) static float dstep; static float complex cnt_last = -I; static float complex dat_last = -I; - static float *cnt_amp; - static float *dat_amp; static float complex *cnt_q; static float complex *dat_q; static struct ddc *cnt_ddc; static struct ddc *dat_ddc; - static struct delay *cnt_delay; - static struct delay *dat_delay; - static short *buff; + static struct buffer *buffer; + static int cnt_delay; + static int dat_delay; + static short *pcm_buff; static int init = 0; if (!init) { init = 1; rate = rate_pcm(pcm); channels = channels_pcm(pcm); - const float step = 1.0 / rate; - #if DN && UP // 320 / 0.088 = 160 / 0.044 = 40000 / 11 = 3636.(36)~ pixels per second for Y, U and V factor_L = 40000; @@ -331,18 +328,19 @@ int demodulate(struct pcm *pcm, float *cnt_freq, float *dat_freq, float *drate) *drate = rate * (float)factor_L / (float)factor_M; dstep = 1.0 / *drate; fprintf(stderr, "using factor of %ld/%ld, working at %.2fhz\n", factor_L, factor_M, *drate); - cnt_amp = malloc(sizeof(float) * factor_M); - dat_amp = malloc(sizeof(float) * factor_M); cnt_q = malloc(sizeof(float complex) * factor_L); dat_q = malloc(sizeof(float complex) * factor_L); // same factor to keep life simple and have accurate horizontal sync - cnt_ddc = alloc_ddc(1200.0, 200.0, step, cnt_taps, factor_L, factor_M, kaiser, 2.0); - dat_ddc = alloc_ddc(1900.0, 800.0, step, dat_taps, factor_L, factor_M, kaiser, 2.0); + cnt_ddc = alloc_ddc(factor_L, factor_M, 1200.0, 200.0, rate, cnt_taps, kaiser, 2.0); + dat_ddc = alloc_ddc(factor_L, factor_M, 1900.0, 800.0, rate, dat_taps, kaiser, 2.0); // delay input by phase shift of other filter to synchronize outputs - cnt_delay = alloc_delay((dat_taps - 1) / (2 * factor_L)); - dat_delay = alloc_delay((cnt_taps - 1) / (2 * factor_L)); + cnt_delay = (dat_taps - 1) / (2 * factor_L); + dat_delay = (cnt_taps - 1) / (2 * factor_L); - buff = (short *)malloc(sizeof(short) * channels * factor_M); + pcm_buff = (short *)malloc(sizeof(short) * channels * factor_M); + + // 0.1 second history + enough room for delay and taps + buffer = alloc_buffer(0.1 * rate + 2 * fmaxf(cnt_delay, dat_delay)); // start immediately below out = factor_L; @@ -350,22 +348,20 @@ int demodulate(struct pcm *pcm, float *cnt_freq, float *dat_freq, float *drate) if (out >= factor_L) { out = 0; - if (!read_pcm(pcm, buff, factor_M)) { + if (!read_pcm(pcm, pcm_buff, factor_M)) { init = 0; - free(buff); + free(pcm_buff); free_ddc(cnt_ddc); free_ddc(dat_ddc); - free(cnt_amp); - free(dat_amp); + free_buffer(buffer); return 0; } - for (int j = 0; j < factor_M; j++) { - float amp = (float)buff[j * channels] / 32767.0; - cnt_amp[j] = do_delay(cnt_delay, amp); - dat_amp[j] = do_delay(dat_delay, amp); - } - do_ddc(cnt_ddc, cnt_amp, cnt_q); - do_ddc(dat_ddc, dat_amp, dat_q); + float *buff = 0; + for (int j = 0; j < factor_M; j++) + buff = do_buffer(buffer, (float)pcm_buff[j * channels] / 32767.0); + + do_ddc(cnt_ddc, buff + cnt_delay, cnt_q); + do_ddc(dat_ddc, buff + dat_delay, dat_q); } *cnt_freq = fclampf(1200.0 + cargf(cnt_q[out] * conjf(cnt_last)) / (2.0 * M_PI * dstep), 1100.0, 1300.0);